home *** CD-ROM | disk | FTP | other *** search
/ Languguage OS 2 / Languguage OS II Version 10-94 (Knowledge Media)(1994).ISO / gnu / calc202a.lha / calc-2.02a / calc-alg-2.el < prev    next >
Lisp/Scheme  |  1993-06-01  |  111KB  |  3,508 lines
  1. ;; Calculator for GNU Emacs, part II [calc-alg-2.el]
  2. ;; Copyright (C) 1990, 1991, 1992, 1993 Free Software Foundation, Inc.
  3. ;; Written by Dave Gillespie, daveg@synaptics.com.
  4.  
  5. ;; This file is part of GNU Emacs.
  6.  
  7. ;; GNU Emacs is distributed in the hope that it will be useful,
  8. ;; but WITHOUT ANY WARRANTY.  No author or distributor
  9. ;; accepts responsibility to anyone for the consequences of using it
  10. ;; or for whether it serves any particular purpose or works at all,
  11. ;; unless he says so in writing.  Refer to the GNU Emacs General Public
  12. ;; License for full details.
  13.  
  14. ;; Everyone is granted permission to copy, modify and redistribute
  15. ;; GNU Emacs, but only under the conditions described in the
  16. ;; GNU Emacs General Public License.   A copy of this license is
  17. ;; supposed to have been given to you along with GNU Emacs so you
  18. ;; can know your rights and responsibilities.  It should be in a
  19. ;; file named COPYING.  Among other things, the copyright notice
  20. ;; and this notice must be preserved on all copies.
  21.  
  22.  
  23.  
  24. ;; This file is autoloaded from calc-ext.el.
  25. (require 'calc-ext)
  26.  
  27. (require 'calc-macs)
  28.  
  29. (defun calc-Need-calc-alg-2 () nil)
  30.  
  31.  
  32. (defun calc-derivative (var num)
  33.   (interactive "sDifferentiate with respect to: \np")
  34.   (calc-slow-wrapper
  35.    (and (< num 0) (error "Order of derivative must be positive"))
  36.    (let ((func (if (calc-is-hyperbolic) 'calcFunc-tderiv 'calcFunc-deriv))
  37.      n expr)
  38.      (if (or (equal var "") (equal var "$"))
  39.      (setq n 2
  40.            expr (calc-top-n 2)
  41.            var (calc-top-n 1))
  42.        (setq var (math-read-expr var))
  43.        (if (eq (car-safe var) 'error)
  44.        (error "Bad format in expression: %s" (nth 1 var)))
  45.        (setq n 1
  46.          expr (calc-top-n 1)))
  47.      (while (>= (setq num (1- num)) 0)
  48.        (setq expr (list func expr var)))
  49.      (calc-enter-result n "derv" expr)))
  50. )
  51.  
  52. (defun calc-integral (var)
  53.   (interactive "sIntegration variable: ")
  54.   (calc-slow-wrapper
  55.    (if (or (equal var "") (equal var "$"))
  56.        (calc-enter-result 2 "intg" (list 'calcFunc-integ
  57.                      (calc-top-n 2)
  58.                      (calc-top-n 1)))
  59.      (let ((var (math-read-expr var)))
  60.        (if (eq (car-safe var) 'error)
  61.        (error "Bad format in expression: %s" (nth 1 var)))
  62.        (calc-enter-result 1 "intg" (list 'calcFunc-integ
  63.                      (calc-top-n 1)
  64.                      var)))))
  65. )
  66.  
  67. (defun calc-num-integral (&optional varname lowname highname)
  68.   (interactive "sIntegration variable: ")
  69.   (calc-tabular-command 'calcFunc-ninteg "Integration" "nint"
  70.             nil varname lowname highname)
  71. )
  72.  
  73. (defun calc-summation (arg &optional varname lowname highname)
  74.   (interactive "P\nsSummation variable: ")
  75.   (calc-tabular-command 'calcFunc-sum "Summation" "sum"
  76.             arg varname lowname highname)
  77. )
  78.  
  79. (defun calc-alt-summation (arg &optional varname lowname highname)
  80.   (interactive "P\nsSummation variable: ")
  81.   (calc-tabular-command 'calcFunc-asum "Summation" "asum"
  82.             arg varname lowname highname)
  83. )
  84.  
  85. (defun calc-product (arg &optional varname lowname highname)
  86.   (interactive "P\nsIndex variable: ")
  87.   (calc-tabular-command 'calcFunc-prod "Index" "prod"
  88.             arg varname lowname highname)
  89. )
  90.  
  91. (defun calc-tabulate (arg &optional varname lowname highname)
  92.   (interactive "P\nsIndex variable: ")
  93.   (calc-tabular-command 'calcFunc-table "Index" "tabl"
  94.             arg varname lowname highname)
  95. )
  96.  
  97. (defun calc-tabular-command (func prompt prefix arg varname lowname highname)
  98.   (calc-slow-wrapper
  99.    (let (var (low nil) (high nil) (step nil) stepname stepnum (num 1) expr)
  100.      (if (consp arg)
  101.      (setq stepnum 1)
  102.        (setq stepnum 0))
  103.      (if (or (equal varname "") (equal varname "$") (null varname))
  104.      (setq high (calc-top-n (+ stepnum 1))
  105.            low (calc-top-n (+ stepnum 2))
  106.            var (calc-top-n (+ stepnum 3))
  107.            num (+ stepnum 4))
  108.        (setq var (if (stringp varname) (math-read-expr varname) varname))
  109.        (if (eq (car-safe var) 'error)
  110.        (error "Bad format in expression: %s" (nth 1 var)))
  111.        (or lowname
  112.        (setq lowname (read-string (concat prompt " variable: " varname
  113.                           ", from: "))))
  114.        (if (or (equal lowname "") (equal lowname "$"))
  115.        (setq high (calc-top-n (+ stepnum 1))
  116.          low (calc-top-n (+ stepnum 2))
  117.          num (+ stepnum 3))
  118.      (setq low (if (stringp lowname) (math-read-expr lowname) lowname))
  119.      (if (eq (car-safe low) 'error)
  120.          (error "Bad format in expression: %s" (nth 1 low)))
  121.      (or highname
  122.          (setq highname (read-string (concat prompt " variable: " varname
  123.                          ", from: " lowname
  124.                          ", to: "))))
  125.      (if (or (equal highname "") (equal highname "$"))
  126.          (setq high (calc-top-n (+ stepnum 1))
  127.            num (+ stepnum 2))
  128.        (setq high (if (stringp highname) (math-read-expr highname)
  129.             highname))
  130.        (if (eq (car-safe high) 'error)
  131.            (error "Bad format in expression: %s" (nth 1 high)))
  132.        (if (consp arg)
  133.            (progn
  134.          (setq stepname (read-string (concat prompt " variable: "
  135.                              varname
  136.                              ", from: " lowname
  137.                              ", to: " highname
  138.                              ", step: ")))
  139.          (if (or (equal stepname "") (equal stepname "$"))
  140.              (setq step (calc-top-n 1)
  141.                num 2)
  142.            (setq step (math-read-expr stepname))
  143.            (if (eq (car-safe step) 'error)
  144.                (error "Bad format in expression: %s"
  145.                   (nth 1 step)))))))))
  146.      (or step
  147.      (if (consp arg)
  148.          (setq step (calc-top-n 1))
  149.        (if arg
  150.            (setq step (prefix-numeric-value arg)))))
  151.      (setq expr (calc-top-n num))
  152.      (calc-enter-result num prefix (append (list func expr var low high)
  153.                        (and step (list step))))))
  154. )
  155.  
  156. (defun calc-solve-for (var)
  157.   (interactive "sVariable to solve for: ")
  158.   (calc-slow-wrapper
  159.    (let ((func (if (calc-is-inverse)
  160.            (if (calc-is-hyperbolic) 'calcFunc-ffinv 'calcFunc-finv)
  161.          (if (calc-is-hyperbolic) 'calcFunc-fsolve 'calcFunc-solve))))
  162.      (if (or (equal var "") (equal var "$"))
  163.      (calc-enter-result 2 "solv" (list func
  164.                        (calc-top-n 2)
  165.                        (calc-top-n 1)))
  166.        (let ((var (if (and (string-match ",\\|[^ ] +[^ ]" var)
  167.                (not (string-match "\\[" var)))
  168.               (math-read-expr (concat "[" var "]"))
  169.             (math-read-expr var))))
  170.      (if (eq (car-safe var) 'error)
  171.          (error "Bad format in expression: %s" (nth 1 var)))
  172.      (calc-enter-result 1 "solv" (list func
  173.                        (calc-top-n 1)
  174.                        var))))))
  175. )
  176.  
  177. (defun calc-poly-roots (var)
  178.   (interactive "sVariable to solve for: ")
  179.   (calc-slow-wrapper
  180.    (if (or (equal var "") (equal var "$"))
  181.        (calc-enter-result 2 "prts" (list 'calcFunc-roots
  182.                      (calc-top-n 2)
  183.                      (calc-top-n 1)))
  184.      (let ((var (if (and (string-match ",\\|[^ ] +[^ ]" var)
  185.              (not (string-match "\\[" var)))
  186.             (math-read-expr (concat "[" var "]"))
  187.           (math-read-expr var))))
  188.        (if (eq (car-safe var) 'error)
  189.        (error "Bad format in expression: %s" (nth 1 var)))
  190.        (calc-enter-result 1 "prts" (list 'calcFunc-roots
  191.                      (calc-top-n 1)
  192.                      var)))))
  193. )
  194.  
  195. (defun calc-taylor (var nterms)
  196.   (interactive "sTaylor expansion variable: \nNNumber of terms: ")
  197.   (calc-slow-wrapper
  198.    (let ((var (math-read-expr var)))
  199.      (if (eq (car-safe var) 'error)
  200.      (error "Bad format in expression: %s" (nth 1 var)))
  201.      (calc-enter-result 1 "tylr" (list 'calcFunc-taylor
  202.                        (calc-top-n 1)
  203.                        var
  204.                        (prefix-numeric-value nterms)))))
  205. )
  206.  
  207.  
  208. (defun math-derivative (expr)   ; uses global values: deriv-var, deriv-total.
  209.   (cond ((equal expr deriv-var)
  210.      1)
  211.     ((or (Math-scalarp expr)
  212.          (eq (car expr) 'sdev)
  213.          (and (eq (car expr) 'var)
  214.           (or (not deriv-total)
  215.               (math-const-var expr)
  216.               (progn
  217.             (math-setup-declarations)
  218.             (memq 'const (nth 1 (or (assq (nth 2 expr)
  219.                               math-decls-cache)
  220.                         math-decls-all)))))))
  221.      0)
  222.     ((eq (car expr) '+)
  223.      (math-add (math-derivative (nth 1 expr))
  224.            (math-derivative (nth 2 expr))))
  225.     ((eq (car expr) '-)
  226.      (math-sub (math-derivative (nth 1 expr))
  227.            (math-derivative (nth 2 expr))))
  228.     ((memq (car expr) '(calcFunc-eq calcFunc-neq calcFunc-lt
  229.                     calcFunc-gt calcFunc-leq calcFunc-geq))
  230.      (list (car expr)
  231.            (math-derivative (nth 1 expr))
  232.            (math-derivative (nth 2 expr))))
  233.     ((eq (car expr) 'neg)
  234.      (math-neg (math-derivative (nth 1 expr))))
  235.     ((eq (car expr) '*)
  236.      (math-add (math-mul (nth 2 expr)
  237.                  (math-derivative (nth 1 expr)))
  238.            (math-mul (nth 1 expr)
  239.                  (math-derivative (nth 2 expr)))))
  240.     ((eq (car expr) '/)
  241.      (math-sub (math-div (math-derivative (nth 1 expr))
  242.                  (nth 2 expr))
  243.            (math-div (math-mul (nth 1 expr)
  244.                        (math-derivative (nth 2 expr)))
  245.                  (math-sqr (nth 2 expr)))))
  246.     ((eq (car expr) '^)
  247.      (let ((du (math-derivative (nth 1 expr)))
  248.            (dv (math-derivative (nth 2 expr))))
  249.        (or (Math-zerop du)
  250.            (setq du (math-mul (nth 2 expr)
  251.                   (math-mul (math-normalize
  252.                          (list '^
  253.                            (nth 1 expr)
  254.                            (math-add (nth 2 expr) -1)))
  255.                         du))))
  256.        (or (Math-zerop dv)
  257.            (setq dv (math-mul (math-normalize
  258.                    (list 'calcFunc-ln (nth 1 expr)))
  259.                   (math-mul expr dv))))
  260.        (math-add du dv)))
  261.     ((eq (car expr) '%)
  262.      (math-derivative (nth 1 expr)))   ; a reasonable definition
  263.     ((eq (car expr) 'vec)
  264.      (math-map-vec 'math-derivative expr))
  265.     ((and (memq (car expr) '(calcFunc-conj calcFunc-re calcFunc-im))
  266.           (= (length expr) 2))
  267.      (list (car expr) (math-derivative (nth 1 expr))))
  268.     ((and (memq (car expr) '(calcFunc-subscr calcFunc-mrow calcFunc-mcol))
  269.           (= (length expr) 3))
  270.      (let ((d (math-derivative (nth 1 expr))))
  271.        (if (math-numberp d)
  272.            0    ; assume x and x_1 are independent vars
  273.          (list (car expr) d (nth 2 expr)))))
  274.     (t (or (and (symbolp (car expr))
  275.             (if (= (length expr) 2)
  276.             (let ((handler (get (car expr) 'math-derivative)))
  277.               (and handler
  278.                    (let ((deriv (math-derivative (nth 1 expr))))
  279.                  (if (Math-zerop deriv)
  280.                      deriv
  281.                    (math-mul (funcall handler (nth 1 expr))
  282.                          deriv)))))
  283.               (let ((handler (get (car expr) 'math-derivative-n)))
  284.             (and handler
  285.                  (funcall handler expr)))))
  286.            (and (not (eq deriv-symb 'pre-expand))
  287.             (let ((exp (math-expand-formula expr)))
  288.               (and exp
  289.                (or (let ((deriv-symb 'pre-expand))
  290.                  (catch 'math-deriv (math-derivative expr)))
  291.                    (math-derivative exp)))))
  292.            (if (or (Math-objvecp expr)
  293.                (eq (car expr) 'var)
  294.                (not (symbolp (car expr))))
  295.            (if deriv-symb
  296.                (throw 'math-deriv nil)
  297.              (list (if deriv-total 'calcFunc-tderiv 'calcFunc-deriv)
  298.                expr
  299.                deriv-var))
  300.          (let ((accum 0)
  301.                (arg expr)
  302.                (n 1)
  303.                derv)
  304.            (while (setq arg (cdr arg))
  305.              (or (Math-zerop (setq derv (math-derivative (car arg))))
  306.              (let ((func (intern (concat (symbol-name (car expr))
  307.                              "'"
  308.                              (if (> n 1)
  309.                              (int-to-string n)
  310.                                ""))))
  311.                    (prop (cond ((= (length expr) 2)
  312.                         'math-derivative-1)
  313.                        ((= (length expr) 3)
  314.                         'math-derivative-2)
  315.                        ((= (length expr) 4)
  316.                         'math-derivative-3)
  317.                        ((= (length expr) 5)
  318.                         'math-derivative-4)
  319.                        ((= (length expr) 6)
  320.                         'math-derivative-5))))
  321.                (setq accum
  322.                  (math-add
  323.                   accum
  324.                   (math-mul
  325.                    derv
  326.                    (let ((handler (get func prop)))
  327.                      (or (and prop handler
  328.                           (apply handler (cdr expr)))
  329.                      (if (and deriv-symb
  330.                           (not (get func
  331.                                 'calc-user-defn)))
  332.                          (throw 'math-deriv nil)
  333.                        (cons func (cdr expr))))))))))
  334.              (setq n (1+ n)))
  335.            accum)))))
  336. )
  337.  
  338. (defun calcFunc-deriv (expr deriv-var &optional deriv-value deriv-symb)
  339.   (let* ((deriv-total nil)
  340.      (res (catch 'math-deriv (math-derivative expr))))
  341.     (or (eq (car-safe res) 'calcFunc-deriv)
  342.     (null res)
  343.     (setq res (math-normalize res)))
  344.     (and res
  345.      (if deriv-value
  346.          (math-expr-subst res deriv-var deriv-value)
  347.        res)))
  348. )
  349.  
  350. (defun calcFunc-tderiv (expr deriv-var &optional deriv-value deriv-symb)
  351.   (math-setup-declarations)
  352.   (let* ((deriv-total t)
  353.      (res (catch 'math-deriv (math-derivative expr))))
  354.     (or (eq (car-safe res) 'calcFunc-tderiv)
  355.     (null res)
  356.     (setq res (math-normalize res)))
  357.     (and res
  358.      (if deriv-value
  359.          (math-expr-subst res deriv-var deriv-value)
  360.        res)))
  361. )
  362.  
  363. (put 'calcFunc-inv\' 'math-derivative-1
  364.      (function (lambda (u) (math-neg (math-div 1 (math-sqr u))))))
  365.  
  366. (put 'calcFunc-sqrt\' 'math-derivative-1
  367.      (function (lambda (u) (math-div 1 (math-mul 2 (list 'calcFunc-sqrt u))))))
  368.  
  369. (put 'calcFunc-deg\' 'math-derivative-1
  370.      (function (lambda (u) (math-div-float '(float 18 1) (math-pi)))))
  371.  
  372. (put 'calcFunc-rad\' 'math-derivative-1
  373.      (function (lambda (u) (math-pi-over-180))))
  374.  
  375. (put 'calcFunc-ln\' 'math-derivative-1
  376.      (function (lambda (u) (math-div 1 u))))
  377.  
  378. (put 'calcFunc-log10\' 'math-derivative-1
  379.      (function (lambda (u)
  380.          (math-div (math-div 1 (math-normalize '(calcFunc-ln 10)))
  381.                u))))
  382.  
  383. (put 'calcFunc-lnp1\' 'math-derivative-1
  384.      (function (lambda (u) (math-div 1 (math-add u 1)))))
  385.  
  386. (put 'calcFunc-log\' 'math-derivative-2
  387.      (function (lambda (x b)
  388.          (and (not (Math-zerop b))
  389.               (let ((lnv (math-normalize
  390.                   (list 'calcFunc-ln b))))
  391.             (math-div 1 (math-mul lnv x)))))))
  392.  
  393. (put 'calcFunc-log\'2 'math-derivative-2
  394.      (function (lambda (x b)
  395.          (let ((lnv (list 'calcFunc-ln b)))
  396.            (math-neg (math-div (list 'calcFunc-log x b)
  397.                        (math-mul lnv b)))))))
  398.  
  399. (put 'calcFunc-exp\' 'math-derivative-1
  400.      (function (lambda (u) (math-normalize (list 'calcFunc-exp u)))))
  401.  
  402. (put 'calcFunc-expm1\' 'math-derivative-1
  403.      (function (lambda (u) (math-normalize (list 'calcFunc-expm1 u)))))
  404.  
  405. (put 'calcFunc-sin\' 'math-derivative-1
  406.      (function (lambda (u) (math-to-radians-2 (math-normalize
  407.                            (list 'calcFunc-cos u))))))
  408.  
  409. (put 'calcFunc-cos\' 'math-derivative-1
  410.      (function (lambda (u) (math-neg (math-to-radians-2
  411.                       (math-normalize
  412.                        (list 'calcFunc-sin u)))))))
  413.  
  414. (put 'calcFunc-tan\' 'math-derivative-1
  415.      (function (lambda (u) (math-to-radians-2
  416.                 (math-div 1 (math-sqr
  417.                      (math-normalize
  418.                       (list 'calcFunc-cos u))))))))
  419.  
  420. (put 'calcFunc-arcsin\' 'math-derivative-1
  421.      (function (lambda (u)
  422.          (math-from-radians-2
  423.           (math-div 1 (math-normalize
  424.                    (list 'calcFunc-sqrt
  425.                      (math-sub 1 (math-sqr u)))))))))
  426.  
  427. (put 'calcFunc-arccos\' 'math-derivative-1
  428.      (function (lambda (u)
  429.          (math-from-radians-2
  430.           (math-div -1 (math-normalize
  431.                 (list 'calcFunc-sqrt
  432.                       (math-sub 1 (math-sqr u)))))))))
  433.  
  434. (put 'calcFunc-arctan\' 'math-derivative-1
  435.      (function (lambda (u) (math-from-radians-2
  436.                 (math-div 1 (math-add 1 (math-sqr u)))))))
  437.  
  438. (put 'calcFunc-sinh\' 'math-derivative-1
  439.      (function (lambda (u) (math-normalize (list 'calcFunc-cosh u)))))
  440.  
  441. (put 'calcFunc-cosh\' 'math-derivative-1
  442.      (function (lambda (u) (math-normalize (list 'calcFunc-sinh u)))))
  443.  
  444. (put 'calcFunc-tanh\' 'math-derivative-1
  445.      (function (lambda (u) (math-div 1 (math-sqr
  446.                     (math-normalize
  447.                      (list 'calcFunc-cosh u)))))))
  448.  
  449. (put 'calcFunc-arcsinh\' 'math-derivative-1
  450.      (function (lambda (u)
  451.          (math-div 1 (math-normalize
  452.                   (list 'calcFunc-sqrt
  453.                     (math-add (math-sqr u) 1)))))))
  454.  
  455. (put 'calcFunc-arccosh\' 'math-derivative-1
  456.      (function (lambda (u)
  457.           (math-div 1 (math-normalize
  458.                    (list 'calcFunc-sqrt
  459.                      (math-add (math-sqr u) -1)))))))
  460.  
  461. (put 'calcFunc-arctanh\' 'math-derivative-1
  462.      (function (lambda (u) (math-div 1 (math-sub 1 (math-sqr u))))))
  463.  
  464. (put 'calcFunc-bern\'2 'math-derivative-2
  465.      (function (lambda (n x)
  466.          (math-mul n (list 'calcFunc-bern (math-add n -1) x)))))
  467.  
  468. (put 'calcFunc-euler\'2 'math-derivative-2
  469.      (function (lambda (n x)
  470.          (math-mul n (list 'calcFunc-euler (math-add n -1) x)))))
  471.  
  472. (put 'calcFunc-gammag\'2 'math-derivative-2
  473.      (function (lambda (a x) (math-deriv-gamma a x 1))))
  474.  
  475. (put 'calcFunc-gammaG\'2 'math-derivative-2
  476.      (function (lambda (a x) (math-deriv-gamma a x -1))))
  477.  
  478. (put 'calcFunc-gammaP\'2 'math-derivative-2
  479.      (function (lambda (a x) (math-deriv-gamma a x
  480.                            (math-div
  481.                         1 (math-normalize
  482.                            (list 'calcFunc-gamma
  483.                              a)))))))
  484.  
  485. (put 'calcFunc-gammaQ\'2 'math-derivative-2
  486.      (function (lambda (a x) (math-deriv-gamma a x
  487.                            (math-div
  488.                         -1 (math-normalize
  489.                             (list 'calcFunc-gamma
  490.                               a)))))))
  491.  
  492. (defun math-deriv-gamma (a x scale)
  493.   (math-mul scale
  494.         (math-mul (math-pow x (math-add a -1))
  495.               (list 'calcFunc-exp (math-neg x))))
  496. )
  497.  
  498. (put 'calcFunc-betaB\' 'math-derivative-3
  499.      (function (lambda (x a b) (math-deriv-beta x a b 1))))
  500.  
  501. (put 'calcFunc-betaI\' 'math-derivative-3
  502.      (function (lambda (x a b) (math-deriv-beta x a b
  503.                         (math-div
  504.                          1 (list 'calcFunc-beta
  505.                              a b))))))
  506.  
  507. (defun math-deriv-beta (x a b scale)
  508.   (math-mul (math-mul (math-pow x (math-add a -1))
  509.               (math-pow (math-sub 1 x) (math-add b -1)))
  510.         scale)
  511. )
  512.  
  513. (put 'calcFunc-erf\' 'math-derivative-1
  514.      (function (lambda (x) (math-div 2
  515.                      (math-mul (list 'calcFunc-exp
  516.                              (math-sqr x))
  517.                            (if calc-symbolic-mode
  518.                            '(calcFunc-sqrt
  519.                              (var pi var-pi))
  520.                          (math-sqrt-pi)))))))
  521.  
  522. (put 'calcFunc-erfc\' 'math-derivative-1
  523.      (function (lambda (x) (math-div -2
  524.                      (math-mul (list 'calcFunc-exp
  525.                              (math-sqr x))
  526.                            (if calc-symbolic-mode
  527.                            '(calcFunc-sqrt
  528.                              (var pi var-pi))
  529.                          (math-sqrt-pi)))))))
  530.  
  531. (put 'calcFunc-besJ\'2 'math-derivative-2
  532.      (function (lambda (v z) (math-div (math-sub (list 'calcFunc-besJ
  533.                                (math-add v -1)
  534.                                z)
  535.                          (list 'calcFunc-besJ
  536.                                (math-add v 1)
  537.                                z))
  538.                        2))))
  539.  
  540. (put 'calcFunc-besY\'2 'math-derivative-2
  541.      (function (lambda (v z) (math-div (math-sub (list 'calcFunc-besY
  542.                                (math-add v -1)
  543.                                z)
  544.                          (list 'calcFunc-besY
  545.                                (math-add v 1)
  546.                                z))
  547.                        2))))
  548.  
  549. (put 'calcFunc-sum 'math-derivative-n
  550.      (function
  551.       (lambda (expr)
  552.     (if (math-expr-contains (cons 'vec (cdr (cdr expr))) deriv-var)
  553.         (throw 'math-deriv nil)
  554.       (cons 'calcFunc-sum
  555.         (cons (math-derivative (nth 1 expr))
  556.               (cdr (cdr expr))))))))
  557.  
  558. (put 'calcFunc-prod 'math-derivative-n
  559.      (function
  560.       (lambda (expr)
  561.     (if (math-expr-contains (cons 'vec (cdr (cdr expr))) deriv-var)
  562.         (throw 'math-deriv nil)
  563.       (math-mul expr
  564.             (cons 'calcFunc-sum
  565.               (cons (math-div (math-derivative (nth 1 expr))
  566.                       (nth 1 expr))
  567.                 (cdr (cdr expr)))))))))
  568.  
  569. (put 'calcFunc-integ 'math-derivative-n
  570.      (function
  571.       (lambda (expr)
  572.     (if (= (length expr) 3)
  573.         (if (equal (nth 2 expr) deriv-var)
  574.         (nth 1 expr)
  575.           (math-normalize
  576.            (list 'calcFunc-integ
  577.              (math-derivative (nth 1 expr))
  578.              (nth 2 expr))))
  579.       (if (= (length expr) 5)
  580.           (let ((lower (math-expr-subst (nth 1 expr) (nth 2 expr)
  581.                         (nth 3 expr)))
  582.             (upper (math-expr-subst (nth 1 expr) (nth 2 expr)
  583.                         (nth 4 expr))))
  584.         (math-add (math-sub (math-mul upper
  585.                           (math-derivative (nth 4 expr)))
  586.                     (math-mul lower
  587.                           (math-derivative (nth 3 expr))))
  588.               (if (equal (nth 2 expr) deriv-var)
  589.                   0
  590.                 (math-normalize
  591.                  (list 'calcFunc-integ
  592.                    (math-derivative (nth 1 expr)) (nth 2 expr)
  593.                    (nth 3 expr) (nth 4 expr)))))))))))
  594.  
  595. (put 'calcFunc-if 'math-derivative-n
  596.      (function
  597.       (lambda (expr)
  598.     (and (= (length expr) 4)
  599.          (list 'calcFunc-if (nth 1 expr)
  600.            (math-derivative (nth 2 expr))
  601.            (math-derivative (nth 3 expr)))))))
  602.  
  603. (put 'calcFunc-subscr 'math-derivative-n
  604.      (function
  605.       (lambda (expr)
  606.     (and (= (length expr) 3)
  607.          (list 'calcFunc-subscr (nth 1 expr)
  608.            (math-derivative (nth 2 expr)))))))
  609.  
  610.  
  611.  
  612.  
  613.  
  614. (setq math-integ-var '(var X ---))
  615. (setq math-integ-var-2 '(var Y ---))
  616. (setq math-integ-vars (list 'f math-integ-var math-integ-var-2))
  617. (setq math-integ-var-list (list math-integ-var))
  618. (setq math-integ-var-list-list (list math-integ-var-list))
  619.  
  620. (defmacro math-tracing-integral (&rest parts)
  621.   (list 'and
  622.     'trace-buffer
  623.     (list 'save-excursion
  624.           '(set-buffer trace-buffer)
  625.           '(goto-char (point-max))
  626.           (list 'and
  627.             '(bolp)
  628.             '(insert (make-string (- math-integral-limit
  629.                          math-integ-level) 32)
  630.                  (format "%2d " math-integ-depth)
  631.                  (make-string math-integ-level 32)))
  632.           ;;(list 'condition-case 'err
  633.             (cons 'insert parts)
  634.         ;;    '(error (insert (prin1-to-string err))))
  635.           '(sit-for 0)))
  636. )
  637.  
  638. ;;; The following wrapper caches results and avoids infinite recursion.
  639. ;;; Each cache entry is: ( A B )          Integral of A is B;
  640. ;;;             ( A N )          Integral of A failed at level N;
  641. ;;;             ( A busy )      Currently working on integral of A;
  642. ;;;             ( A parts )      Currently working, integ-by-parts;
  643. ;;;             ( A parts2 )      Currently working, integ-by-parts;
  644. ;;;             ( A cancelled )  Ignore this cache entry;
  645. ;;;             ( A [B] )        Same result as for cur-record = B.
  646. (defun math-integral (expr &optional simplify same-as-above)
  647.   (let* ((simp cur-record)
  648.      (cur-record (assoc expr math-integral-cache))
  649.      (math-integ-depth (1+ math-integ-depth))
  650.      (val 'cancelled))
  651.     (math-tracing-integral "Integrating "
  652.                (math-format-value expr 1000)
  653.                "...\n")
  654.     (and cur-record
  655.      (progn
  656.        (math-tracing-integral "Found "
  657.                   (math-format-value (nth 1 cur-record) 1000))
  658.        (and (consp (nth 1 cur-record))
  659.         (math-replace-integral-parts cur-record))
  660.        (math-tracing-integral " => "
  661.                   (math-format-value (nth 1 cur-record) 1000)
  662.                   "\n")))
  663.     (or (and cur-record
  664.          (not (eq (nth 1 cur-record) 'cancelled))
  665.          (or (not (integerp (nth 1 cur-record)))
  666.          (>= (nth 1 cur-record) math-integ-level)))
  667.     (and (math-integral-contains-parts expr)
  668.          (progn
  669.            (setq val nil)
  670.            t))
  671.     (unwind-protect
  672.         (progn
  673.           (let (math-integ-msg)
  674.         (if (eq calc-display-working-message 'lots)
  675.             (progn
  676.               (calc-set-command-flag 'clear-message)
  677.               (setq math-integ-msg (format
  678.                         "Working... Integrating %s"
  679.                         (math-format-flat-expr expr 0)))
  680.               (message math-integ-msg)))
  681.         (if cur-record
  682.             (setcar (cdr cur-record)
  683.                 (if same-as-above (vector simp) 'busy))
  684.           (setq cur-record
  685.             (list expr (if same-as-above (vector simp) 'busy))
  686.             math-integral-cache (cons cur-record
  687.                           math-integral-cache)))
  688.         (if (eq simplify 'yes)
  689.             (progn
  690.               (math-tracing-integral "Simplifying...")
  691.               (setq simp (math-simplify expr))
  692.               (setq val (if (equal simp expr)
  693.                     (progn
  694.                       (math-tracing-integral " no change\n")
  695.                       (math-do-integral expr))
  696.                   (math-tracing-integral " simplified\n")
  697.                   (math-integral simp 'no t))))
  698.           (or (setq val (math-do-integral expr))
  699.               (eq simplify 'no)
  700.               (let ((simp (math-simplify expr)))
  701.             (or (equal simp expr)
  702.                 (progn
  703.                   (math-tracing-integral "Trying again after "
  704.                              "simplification...\n")
  705.                   (setq val (math-integral simp 'no t))))))))
  706.           (if (eq calc-display-working-message 'lots)
  707.           (message math-integ-msg)))
  708.       (setcar (cdr cur-record) (or val
  709.                        (if (or math-enable-subst
  710.                            (not math-any-substs))
  711.                        math-integ-level
  712.                      'cancelled)))))
  713.     (setq val cur-record)
  714.     (while (vectorp (nth 1 val))
  715.       (setq val (aref (nth 1 val) 0)))
  716.     (setq val (if (memq (nth 1 val) '(parts parts2))
  717.           (progn
  718.             (setcar (cdr val) 'parts2)
  719.             (list 'var 'PARTS val))
  720.         (and (consp (nth 1 val))
  721.              (nth 1 val))))
  722.     (math-tracing-integral "Integral of "
  723.                (math-format-value expr 1000)
  724.                "  is  "
  725.                (math-format-value val 1000)
  726.                "\n")
  727.     val)
  728. )
  729. (defvar math-integral-cache nil)
  730. (defvar math-integral-cache-state nil)
  731.  
  732. (defun math-integral-contains-parts (expr)
  733.   (if (Math-primp expr)
  734.       (and (eq (car-safe expr) 'var)
  735.        (eq (nth 1 expr) 'PARTS)
  736.        (listp (nth 2 expr)))
  737.     (while (and (setq expr (cdr expr))
  738.         (not (math-integral-contains-parts (car expr)))))
  739.     expr)
  740. )
  741.  
  742. (defun math-replace-integral-parts (expr)
  743.   (or (Math-primp expr)
  744.       (while (setq expr (cdr expr))
  745.     (and (consp (car expr))
  746.          (if (eq (car (car expr)) 'var)
  747.          (and (eq (nth 1 (car expr)) 'PARTS)
  748.               (consp (nth 2 (car expr)))
  749.               (if (listp (nth 1 (nth 2 (car expr))))
  750.               (progn
  751.                 (setcar expr (nth 1 (nth 2 (car expr))))
  752.                 (math-replace-integral-parts (cons 'foo expr)))
  753.             (setcar (cdr cur-record) 'cancelled)))
  754.            (math-replace-integral-parts (car expr))))))
  755. )
  756.  
  757. (defun math-do-integral (expr)
  758.   (let (t1 t2)
  759.     (or (cond ((not (math-expr-contains expr math-integ-var))
  760.            (math-mul expr math-integ-var))
  761.           ((equal expr math-integ-var)
  762.            (math-div (math-sqr expr) 2))
  763.           ((eq (car expr) '+)
  764.            (and (setq t1 (math-integral (nth 1 expr)))
  765.             (setq t2 (math-integral (nth 2 expr)))
  766.             (math-add t1 t2)))
  767.           ((eq (car expr) '-)
  768.            (and (setq t1 (math-integral (nth 1 expr)))
  769.             (setq t2 (math-integral (nth 2 expr)))
  770.             (math-sub t1 t2)))
  771.           ((eq (car expr) 'neg)
  772.            (and (setq t1 (math-integral (nth 1 expr)))
  773.             (math-neg t1)))
  774.           ((eq (car expr) '*)
  775.            (cond ((not (math-expr-contains (nth 1 expr) math-integ-var))
  776.               (and (setq t1 (math-integral (nth 2 expr)))
  777.                (math-mul (nth 1 expr) t1)))
  778.              ((not (math-expr-contains (nth 2 expr) math-integ-var))
  779.               (and (setq t1 (math-integral (nth 1 expr)))
  780.                (math-mul t1 (nth 2 expr))))
  781.              ((memq (car-safe (nth 1 expr)) '(+ -))
  782.               (math-integral (list (car (nth 1 expr))
  783.                        (math-mul (nth 1 (nth 1 expr))
  784.                              (nth 2 expr))
  785.                        (math-mul (nth 2 (nth 1 expr))
  786.                              (nth 2 expr)))
  787.                      'yes t))
  788.              ((memq (car-safe (nth 2 expr)) '(+ -))
  789.               (math-integral (list (car (nth 2 expr))
  790.                        (math-mul (nth 1 (nth 2 expr))
  791.                              (nth 1 expr))
  792.                        (math-mul (nth 2 (nth 2 expr))
  793.                              (nth 1 expr)))
  794.                      'yes t))))
  795.           ((eq (car expr) '/)
  796.            (cond ((and (not (math-expr-contains (nth 1 expr)
  797.                             math-integ-var))
  798.                (not (math-equal-int (nth 1 expr) 1)))
  799.               (and (setq t1 (math-integral (math-div 1 (nth 2 expr))))
  800.                (math-mul (nth 1 expr) t1)))
  801.              ((not (math-expr-contains (nth 2 expr) math-integ-var))
  802.               (and (setq t1 (math-integral (nth 1 expr)))
  803.                (math-div t1 (nth 2 expr))))
  804.              ((and (eq (car-safe (nth 1 expr)) '*)
  805.                (not (math-expr-contains (nth 1 (nth 1 expr))
  806.                             math-integ-var)))
  807.               (and (setq t1 (math-integral
  808.                      (math-div (nth 2 (nth 1 expr))
  809.                            (nth 2 expr))))
  810.                (math-mul t1 (nth 1 (nth 1 expr)))))
  811.              ((and (eq (car-safe (nth 1 expr)) '*)
  812.                (not (math-expr-contains (nth 2 (nth 1 expr))
  813.                             math-integ-var)))
  814.               (and (setq t1 (math-integral
  815.                      (math-div (nth 1 (nth 1 expr))
  816.                            (nth 2 expr))))
  817.                (math-mul t1 (nth 2 (nth 1 expr)))))
  818.              ((and (eq (car-safe (nth 2 expr)) '*)
  819.                (not (math-expr-contains (nth 1 (nth 2 expr))
  820.                             math-integ-var)))
  821.               (and (setq t1 (math-integral
  822.                      (math-div (nth 1 expr)
  823.                            (nth 2 (nth 2 expr)))))
  824.                (math-div t1 (nth 1 (nth 2 expr)))))
  825.              ((and (eq (car-safe (nth 2 expr)) '*)
  826.                (not (math-expr-contains (nth 2 (nth 2 expr))
  827.                             math-integ-var)))
  828.               (and (setq t1 (math-integral
  829.                      (math-div (nth 1 expr)
  830.                            (nth 1 (nth 2 expr)))))
  831.                (math-div t1 (nth 2 (nth 2 expr)))))
  832.              ((eq (car-safe (nth 2 expr)) 'calcFunc-exp)
  833.               (math-integral
  834.                (math-mul (nth 1 expr)
  835.                  (list 'calcFunc-exp
  836.                        (math-neg (nth 1 (nth 2 expr)))))))))
  837.           ((eq (car expr) '^)
  838.            (cond ((not (math-expr-contains (nth 1 expr) math-integ-var))
  839.               (or (and (setq t1 (math-is-polynomial (nth 2 expr)
  840.                                 math-integ-var 1))
  841.                    (math-div expr
  842.                      (math-mul (nth 1 t1)
  843.                            (math-normalize
  844.                             (list 'calcFunc-ln
  845.                               (nth 1 expr))))))
  846.               (math-integral
  847.                (list 'calcFunc-exp
  848.                  (math-mul (nth 2 expr)
  849.                        (math-normalize
  850.                         (list 'calcFunc-ln
  851.                           (nth 1 expr)))))
  852.                'yes t)))
  853.              ((not (math-expr-contains (nth 2 expr) math-integ-var))
  854.               (if (and (integerp (nth 2 expr)) (< (nth 2 expr) 0))
  855.               (math-integral
  856.                (list '/ 1 (math-pow (nth 1 expr) (- (nth 2 expr))))
  857.                nil t)
  858.             (or (and (setq t1 (math-is-polynomial (nth 1 expr)
  859.                                   math-integ-var
  860.                                   1))
  861.                  (setq t2 (math-add (nth 2 expr) 1))
  862.                  (math-div (math-pow (nth 1 expr) t2)
  863.                        (math-mul t2 (nth 1 t1))))
  864.                 (and (Math-negp (nth 2 expr))
  865.                  (math-integral
  866.                   (math-div 1
  867.                         (math-pow (nth 1 expr)
  868.                               (math-neg
  869.                                (nth 2 expr))))
  870.                   nil t))
  871.                 nil))))))
  872.  
  873.     ;; Integral of a polynomial.
  874.     (and (setq t1 (math-is-polynomial expr math-integ-var 20))
  875.          (let ((accum 0)
  876.            (n 1))
  877.            (while t1
  878.          (if (setq accum (math-add accum
  879.                        (math-div (math-mul (car t1)
  880.                                    (math-pow
  881.                                 math-integ-var
  882.                                 n))
  883.                              n))
  884.                t1 (cdr t1))
  885.              (setq n (1+ n))))
  886.            accum))
  887.  
  888.     ;; Try looking it up!
  889.     (cond ((= (length expr) 2)
  890.            (and (symbolp (car expr))
  891.             (setq t1 (get (car expr) 'math-integral))
  892.             (progn
  893.               (while (and t1
  894.                   (not (setq t2 (funcall (car t1)
  895.                              (nth 1 expr)))))
  896.             (setq t1 (cdr t1)))
  897.               (and t2 (math-normalize t2)))))
  898.           ((= (length expr) 3)
  899.            (and (symbolp (car expr))
  900.             (setq t1 (get (car expr) 'math-integral-2))
  901.             (progn
  902.               (while (and t1
  903.                   (not (setq t2 (funcall (car t1)
  904.                              (nth 1 expr)
  905.                              (nth 2 expr)))))
  906.             (setq t1 (cdr t1)))
  907.               (and t2 (math-normalize t2))))))
  908.  
  909.     ;; Integral of a rational function.
  910.     (and (math-ratpoly-p expr math-integ-var)
  911.          (setq t1 (calcFunc-apart expr math-integ-var))
  912.          (not (equal t1 expr))
  913.          (math-integral t1))
  914.  
  915.     ;; Try user-defined integration rules.
  916.     (and has-rules
  917.          (let ((math-old-integ (symbol-function 'calcFunc-integ))
  918.            (input (list 'calcFunc-integtry expr math-integ-var))
  919.            res part)
  920.            (unwind-protect
  921.            (progn
  922.              (fset 'calcFunc-integ 'math-sub-integration)
  923.              (setq res (math-rewrite input
  924.                          '(var IntegRules var-IntegRules)
  925.                          1))
  926.              (fset 'calcFunc-integ math-old-integ)
  927.              (and (not (equal res input))
  928.               (if (setq part (math-expr-calls
  929.                       res '(calcFunc-integsubst)))
  930.                   (and (memq (length part) '(3 4 5))
  931.                    (let ((parts (mapcar
  932.                          (function
  933.                           (lambda (x)
  934.                             (math-expr-subst
  935.                              x (nth 2 part)
  936.                              math-integ-var)))
  937.                          (cdr part))))
  938.                      (math-integrate-by-substitution
  939.                       expr (car parts) t
  940.                       (or (nth 2 parts)
  941.                       (list 'calcFunc-integfailed
  942.                         math-integ-var))
  943.                       (nth 3 parts))))
  944.                 (if (not (math-expr-calls res
  945.                               '(calcFunc-integtry
  946.                             calcFunc-integfailed)))
  947.                 res))))
  948.          (fset 'calcFunc-integ math-old-integ))))
  949.  
  950.     ;; See if the function is a symbolic derivative.
  951.     (and (string-match "'" (symbol-name (car expr)))
  952.          (let ((name (symbol-name (car expr)))
  953.            (p expr) (n 0) (which nil) (bad nil))
  954.            (while (setq n (1+ n) p (cdr p))
  955.          (if (equal (car p) math-integ-var)
  956.              (if which (setq bad t) (setq which n))
  957.            (if (math-expr-contains (car p) math-integ-var)
  958.                (setq bad t))))
  959.            (and which (not bad)
  960.             (let ((prime (if (= which 1) "'" (format "'%d" which))))
  961.               (and (string-match (concat prime "\\('['0-9]*\\|$\\)")
  962.                      name)
  963.                (cons (intern
  964.                   (concat
  965.                    (substring name 0 (match-beginning 0))
  966.                    (substring name (+ (match-beginning 0)
  967.                               (length prime)))))
  968.                  (cdr expr)))))))
  969.  
  970.     ;; Try transformation methods (parts, substitutions).
  971.     (and (> math-integ-level 0)
  972.          (math-do-integral-methods expr))
  973.  
  974.     ;; Try expanding the function's definition.
  975.     (let ((res (math-expand-formula expr)))
  976.       (and res
  977.            (math-integral res)))))
  978. )
  979.  
  980. (defun math-sub-integration (expr &rest rest)
  981.   (or (if (or (not rest)
  982.           (and (< math-integ-level math-integral-limit)
  983.            (eq (car rest) math-integ-var)))
  984.       (math-integral expr)
  985.     (let ((res (apply math-old-integ expr rest)))
  986.       (and (or (= math-integ-level math-integral-limit)
  987.            (not (math-expr-calls res 'calcFunc-integ)))
  988.            res)))
  989.       (list 'calcFunc-integfailed expr))
  990. )
  991.  
  992. (defun math-do-integral-methods (expr)
  993.   (let ((so-far math-integ-var-list-list)
  994.     rat-in)
  995.  
  996.     ;; Integration by substitution, for various likely sub-expressions.
  997.     ;; (In first pass, we look only for sub-exprs that are linear in X.)
  998.     (or (if math-enable-subst
  999.         (math-integ-try-substitutions expr)
  1000.       (math-integ-try-linear-substitutions expr))
  1001.  
  1002.     ;; If function has sines and cosines, try tan(x/2) substitution.
  1003.     (and (let ((p (setq rat-in (math-expr-rational-in expr))))
  1004.            (while (and p
  1005.                (memq (car (car p)) '(calcFunc-sin
  1006.                          calcFunc-cos
  1007.                          calcFunc-tan))
  1008.                (equal (nth 1 (car p)) math-integ-var))
  1009.          (setq p (cdr p)))
  1010.            (null p))
  1011.          (or (and (math-integ-parts-easy expr)
  1012.               (math-integ-try-parts expr t))
  1013.          (math-integrate-by-good-substitution
  1014.           expr (list 'calcFunc-tan (math-div math-integ-var 2)))))
  1015.  
  1016.     ;; If function has sinh and cosh, try tanh(x/2) substitution.
  1017.     (and (let ((p rat-in))
  1018.            (while (and p
  1019.                (memq (car (car p)) '(calcFunc-sinh
  1020.                          calcFunc-cosh
  1021.                          calcFunc-tanh
  1022.                          calcFunc-exp))
  1023.                (equal (nth 1 (car p)) math-integ-var))
  1024.          (setq p (cdr p)))
  1025.            (null p))
  1026.          (or (and (math-integ-parts-easy expr)
  1027.               (math-integ-try-parts expr t))
  1028.          (math-integrate-by-good-substitution
  1029.           expr (list 'calcFunc-tanh (math-div math-integ-var 2)))))
  1030.  
  1031.     ;; If function has square roots, try sin, tan, or sec substitution.
  1032.     (and (let ((p rat-in))
  1033.            (setq t1 nil)
  1034.            (while (and p
  1035.                (or (equal (car p) math-integ-var)
  1036.                    (and (eq (car (car p)) 'calcFunc-sqrt)
  1037.                     (setq t1 (math-is-polynomial
  1038.                           (nth 1 (setq t2 (car p)))
  1039.                           math-integ-var 2)))))
  1040.          (setq p (cdr p)))
  1041.            (and (null p) t1))
  1042.          (if (cdr (cdr t1))
  1043.          (if (math-guess-if-neg (nth 2 t1))
  1044.              (let* ((c (math-sqrt (math-neg (nth 2 t1))))
  1045.                 (d (math-div (nth 1 t1) (math-mul -2 c)))
  1046.                 (a (math-sqrt (math-add (car t1) (math-sqr d)))))
  1047.                (math-integrate-by-good-substitution
  1048.             expr (list 'calcFunc-arcsin
  1049.                    (math-div-thru
  1050.                     (math-add (math-mul c math-integ-var) d)
  1051.                     a))))
  1052.            (let* ((c (math-sqrt (nth 2 t1)))
  1053.               (d (math-div (nth 1 t1) (math-mul 2 c)))
  1054.               (aa (math-sub (car t1) (math-sqr d))))
  1055.              (if (and nil (not (and (eq d 0) (eq c 1))))
  1056.              (math-integrate-by-good-substitution
  1057.               expr (math-add (math-mul c math-integ-var) d))
  1058.                (if (math-guess-if-neg aa)
  1059.                (math-integrate-by-good-substitution
  1060.                 expr (list 'calcFunc-arccosh
  1061.                        (math-div-thru
  1062.                     (math-add (math-mul c math-integ-var)
  1063.                           d)
  1064.                     (math-sqrt (math-neg aa)))))
  1065.              (math-integrate-by-good-substitution
  1066.               expr (list 'calcFunc-arcsinh
  1067.                      (math-div-thru
  1068.                       (math-add (math-mul c math-integ-var)
  1069.                         d)
  1070.                       (math-sqrt aa))))))))
  1071.            (math-integrate-by-good-substitution expr t2)) )
  1072.  
  1073.     ;; Try integration by parts.
  1074.     (math-integ-try-parts expr)
  1075.  
  1076.     ;; Give up.
  1077.     nil))
  1078. )
  1079.  
  1080. (defun math-integ-parts-easy (expr)
  1081.   (cond ((Math-primp expr) t)
  1082.     ((memq (car expr) '(+ - *))
  1083.      (and (math-integ-parts-easy (nth 1 expr))
  1084.           (math-integ-parts-easy (nth 2 expr))))
  1085.     ((eq (car expr) '/)
  1086.      (and (math-integ-parts-easy (nth 1 expr))
  1087.           (math-atomic-factorp (nth 2 expr))))
  1088.     ((eq (car expr) '^)
  1089.      (and (natnump (nth 2 expr))
  1090.           (math-integ-parts-easy (nth 1 expr))))
  1091.     ((eq (car expr) 'neg)
  1092.      (math-integ-parts-easy (nth 1 expr)))
  1093.     (t t))
  1094. )
  1095.  
  1096. (defun math-integ-try-parts (expr &optional math-good-parts)
  1097.   ;; Integration by parts:
  1098.   ;;   integ(f(x) g(x),x) = f(x) h(x) - integ(h(x) f'(x),x)
  1099.   ;;     where h(x) = integ(g(x),x).
  1100.   (or (let ((exp (calcFunc-expand expr)))
  1101.     (and (not (equal exp expr))
  1102.          (math-integral exp)))
  1103.       (and (eq (car expr) '*)
  1104.        (let ((first-bad (or (math-polynomial-p (nth 1 expr)
  1105.                            math-integ-var)
  1106.                 (equal (nth 2 expr) math-prev-parts-v))))
  1107.          (or (and first-bad   ; so try this one first
  1108.               (math-integrate-by-parts (nth 1 expr) (nth 2 expr)))
  1109.          (math-integrate-by-parts (nth 2 expr) (nth 1 expr))
  1110.          (and (not first-bad)
  1111.               (math-integrate-by-parts (nth 1 expr) (nth 2 expr))))))
  1112.       (and (eq (car expr) '/)
  1113.        (math-expr-contains (nth 1 expr) math-integ-var)
  1114.        (let ((recip (math-div 1 (nth 2 expr))))
  1115.          (or (math-integrate-by-parts (nth 1 expr) recip)
  1116.          (math-integrate-by-parts recip (nth 1 expr)))))
  1117.       (and (eq (car expr) '^)
  1118.        (math-integrate-by-parts (math-pow (nth 1 expr)
  1119.                           (math-sub (nth 2 expr) 1))
  1120.                     (nth 1 expr))))
  1121. )
  1122.  
  1123. (defun math-integrate-by-parts (u vprime)
  1124.   (let ((math-integ-level (if (or math-good-parts
  1125.                   (math-polynomial-p u math-integ-var))
  1126.                   math-integ-level
  1127.                 (1- math-integ-level)))
  1128.     (math-doing-parts t)
  1129.     v temp)
  1130.     (and (>= math-integ-level 0)
  1131.      (unwind-protect
  1132.          (progn
  1133.            (setcar (cdr cur-record) 'parts)
  1134.            (math-tracing-integral "Integrating by parts, u = "
  1135.                       (math-format-value u 1000)
  1136.                       ", v' = "
  1137.                       (math-format-value vprime 1000)
  1138.                       "\n")
  1139.            (and (setq v (math-integral vprime))
  1140.             (setq temp (calcFunc-deriv u math-integ-var nil t))
  1141.             (setq temp (let ((math-prev-parts-v v))
  1142.                  (math-integral (math-mul v temp) 'yes)))
  1143.             (setq temp (math-sub (math-mul u v) temp))
  1144.             (if (eq (nth 1 cur-record) 'parts)
  1145.             (calcFunc-expand temp)
  1146.               (setq v (list 'var 'PARTS cur-record)
  1147.                 var-thing (list 'vec (math-sub v temp) v)
  1148.                 temp (let (calc-next-why)
  1149.                    (math-solve-for (math-sub v temp) 0 v nil)))
  1150.               (and temp (not (integerp temp))
  1151.                (math-simplify-extended temp)))))
  1152.        (setcar (cdr cur-record) 'busy))))
  1153. )
  1154.  
  1155. ;;; This tries two different formulations, hoping the algebraic simplifier
  1156. ;;; will be strong enough to handle at least one.
  1157. (defun math-integrate-by-substitution (expr u &optional user uinv uinvprime)
  1158.   (and (> math-integ-level 0)
  1159.        (let ((math-integ-level (max (- math-integ-level 2) 0)))
  1160.      (math-integrate-by-good-substitution expr u user uinv uinvprime)))
  1161. )
  1162.  
  1163. (defun math-integrate-by-good-substitution (expr u &optional user
  1164.                          uinv uinvprime)
  1165.   (let ((math-living-dangerously t)
  1166.     deriv temp)
  1167.     (and (setq uinv (if uinv
  1168.             (math-expr-subst uinv math-integ-var
  1169.                      math-integ-var-2)
  1170.               (let (calc-next-why)
  1171.             (math-solve-for u
  1172.                     math-integ-var-2
  1173.                     math-integ-var nil))))
  1174.      (progn
  1175.        (math-tracing-integral "Integrating by substitution, u = "
  1176.                   (math-format-value u 1000)
  1177.                   "\n")
  1178.        (or (and (setq deriv (calcFunc-deriv u
  1179.                         math-integ-var nil
  1180.                         (not user)))
  1181.             (setq temp (math-integral (math-expr-subst
  1182.                            (math-expr-subst
  1183.                         (math-expr-subst
  1184.                          (math-div expr deriv)
  1185.                          u
  1186.                          math-integ-var-2)
  1187.                         math-integ-var
  1188.                         uinv)
  1189.                            math-integ-var-2
  1190.                            math-integ-var)
  1191.                           'yes)))
  1192.            (and (setq deriv (or uinvprime
  1193.                     (calcFunc-deriv uinv
  1194.                             math-integ-var-2
  1195.                             math-integ-var
  1196.                             (not user))))
  1197.             (setq temp (math-integral (math-mul
  1198.                            (math-expr-subst
  1199.                         (math-expr-subst
  1200.                          (math-expr-subst
  1201.                           expr
  1202.                           u
  1203.                           math-integ-var-2)
  1204.                          math-integ-var
  1205.                          uinv)
  1206.                         math-integ-var-2
  1207.                         math-integ-var)
  1208.                            deriv)
  1209.                           'yes)))))
  1210.      (math-simplify-extended
  1211.       (math-expr-subst temp math-integ-var u))))
  1212. )
  1213.  
  1214. ;;; Look for substitutions of the form u = a x + b.
  1215. (defun math-integ-try-linear-substitutions (sub-expr)
  1216.   (and (not (Math-primp sub-expr))
  1217.        (or (and (not (memq (car sub-expr) '(+ - * / neg)))
  1218.         (not (and (eq (car sub-expr) '^)
  1219.               (integerp (nth 2 sub-expr))))
  1220.         (math-expr-contains sub-expr math-integ-var)
  1221.         (let ((res nil))
  1222.           (while (and (setq sub-expr (cdr sub-expr))
  1223.                   (or (not (math-linear-in (car sub-expr)
  1224.                                math-integ-var))
  1225.                   (assoc (car sub-expr) so-far)
  1226.                   (progn
  1227.                     (setq so-far (cons (list (car sub-expr))
  1228.                                so-far))
  1229.                     (not (setq res
  1230.                            (math-integrate-by-substitution
  1231.                         expr (car sub-expr))))))))
  1232.           res))
  1233.        (let ((res nil))
  1234.          (while (and (setq sub-expr (cdr sub-expr))
  1235.              (not (setq res (math-integ-try-linear-substitutions
  1236.                      (car sub-expr))))))
  1237.          res)))
  1238. )
  1239.  
  1240. ;;; Recursively try different substitutions based on various sub-expressions.
  1241. (defun math-integ-try-substitutions (sub-expr &optional allow-rat)
  1242.   (and (not (Math-primp sub-expr))
  1243.        (not (assoc sub-expr so-far))
  1244.        (math-expr-contains sub-expr math-integ-var)
  1245.        (or (and (if (and (not (memq (car sub-expr) '(+ - * / neg)))
  1246.              (not (and (eq (car sub-expr) '^)
  1247.                    (integerp (nth 2 sub-expr)))))
  1248.             (setq allow-rat t)
  1249.           (prog1 allow-rat (setq allow-rat nil)))
  1250.         (not (eq sub-expr expr))
  1251.         (or (math-integrate-by-substitution expr sub-expr)
  1252.             (and (eq (car sub-expr) '^)
  1253.              (integerp (nth 2 sub-expr))
  1254.              (< (nth 2 sub-expr) 0)
  1255.              (math-integ-try-substitutions
  1256.               (math-pow (nth 1 sub-expr) (- (nth 2 sub-expr)))
  1257.               t))))
  1258.        (let ((res nil))
  1259.          (setq so-far (cons (list sub-expr) so-far))
  1260.          (while (and (setq sub-expr (cdr sub-expr))
  1261.              (not (setq res (math-integ-try-substitutions
  1262.                      (car sub-expr) allow-rat)))))
  1263.          res)))
  1264. )
  1265.  
  1266. (defun math-expr-rational-in (expr)
  1267.   (let ((parts nil))
  1268.     (math-expr-rational-in-rec expr)
  1269.     (mapcar 'car parts))
  1270. )
  1271.  
  1272. (defun math-expr-rational-in-rec (expr)
  1273.   (cond ((Math-primp expr)
  1274.      (and (equal expr math-integ-var)
  1275.           (not (assoc expr parts))
  1276.           (setq parts (cons (list expr) parts))))
  1277.     ((or (memq (car expr) '(+ - * / neg))
  1278.          (and (eq (car expr) '^) (integerp (nth 2 expr))))
  1279.      (math-expr-rational-in-rec (nth 1 expr))
  1280.      (and (nth 2 expr) (math-expr-rational-in-rec (nth 2 expr))))
  1281.     ((and (eq (car expr) '^)
  1282.           (eq (math-quarter-integer (nth 2 expr)) 2))
  1283.      (math-expr-rational-in-rec (list 'calcFunc-sqrt (nth 1 expr))))
  1284.     (t
  1285.      (and (not (assoc expr parts))
  1286.           (math-expr-contains expr math-integ-var)
  1287.           (setq parts (cons (list expr) parts)))))
  1288. )
  1289.  
  1290. (defun math-expr-calls (expr funcs &optional arg-contains)
  1291.   (if (consp expr)
  1292.       (if (or (memq (car expr) funcs)
  1293.           (and (eq (car expr) '^) (eq (car funcs) 'calcFunc-sqrt)
  1294.            (eq (math-quarter-integer (nth 2 expr)) 2)))
  1295.       (and (or (not arg-contains)
  1296.            (math-expr-contains expr arg-contains))
  1297.            expr)
  1298.     (and (not (Math-primp expr))
  1299.          (let ((res nil))
  1300.            (while (and (setq expr (cdr expr))
  1301.                (not (setq res (math-expr-calls
  1302.                        (car expr) funcs arg-contains)))))
  1303.            res))))
  1304. )
  1305.  
  1306. (defun math-fix-const-terms (expr except-vars)
  1307.   (cond ((not (math-expr-depends expr except-vars)) 0)
  1308.     ((Math-primp expr) expr)
  1309.     ((eq (car expr) '+)
  1310.      (math-add (math-fix-const-terms (nth 1 expr) except-vars)
  1311.            (math-fix-const-terms (nth 2 expr) except-vars)))
  1312.     ((eq (car expr) '-)
  1313.      (math-sub (math-fix-const-terms (nth 1 expr) except-vars)
  1314.            (math-fix-const-terms (nth 2 expr) except-vars)))
  1315.     (t expr))
  1316. )
  1317.  
  1318. ;; Command for debugging the Calculator's symbolic integrator.
  1319. (defun calc-dump-integral-cache (&optional arg)
  1320.   (interactive "P")
  1321.   (let ((buf (current-buffer)))
  1322.     (unwind-protect
  1323.     (let ((p math-integral-cache)
  1324.           cur-record)
  1325.       (display-buffer (get-buffer-create "*Integral Cache*")) 
  1326.       (set-buffer (get-buffer "*Integral Cache*"))
  1327.       (erase-buffer)
  1328.       (while p
  1329.         (setq cur-record (car p))
  1330.         (or arg (math-replace-integral-parts cur-record))
  1331.         (insert (math-format-flat-expr (car cur-record) 0)
  1332.             " --> "
  1333.             (if (symbolp (nth 1 cur-record))
  1334.             (concat "(" (symbol-name (nth 1 cur-record)) ")")
  1335.               (math-format-flat-expr (nth 1 cur-record) 0))
  1336.             "\n")
  1337.         (setq p (cdr p)))
  1338.       (goto-char (point-min)))
  1339.       (set-buffer buf)))
  1340. )
  1341.  
  1342. (defun math-try-integral (expr)
  1343.   (let ((math-integ-level math-integral-limit)
  1344.     (math-integ-depth 0)
  1345.     (math-integ-msg "Working...done")
  1346.     (cur-record nil)   ; a technicality
  1347.     (math-integrating t)
  1348.     (calc-prefer-frac t)
  1349.     (calc-symbolic-mode t)
  1350.     (has-rules (calc-has-rules 'var-IntegRules)))
  1351.     (or (math-integral expr 'yes)
  1352.     (and math-any-substs
  1353.          (setq math-enable-subst t)
  1354.          (math-integral expr 'yes))
  1355.     (and (> math-max-integral-limit math-integral-limit)
  1356.          (setq math-integral-limit math-max-integral-limit
  1357.            math-integ-level math-integral-limit)
  1358.          (math-integral expr 'yes))))
  1359. )
  1360.  
  1361. (defun calcFunc-integ (expr var &optional low high)
  1362.   (cond
  1363.    ;; Do these even if the parts turn out not to be integrable.
  1364.    ((eq (car-safe expr) '+)
  1365.     (math-add (calcFunc-integ (nth 1 expr) var low high)
  1366.           (calcFunc-integ (nth 2 expr) var low high)))
  1367.    ((eq (car-safe expr) '-)
  1368.     (math-sub (calcFunc-integ (nth 1 expr) var low high)
  1369.           (calcFunc-integ (nth 2 expr) var low high)))
  1370.    ((eq (car-safe expr) 'neg)
  1371.     (math-neg (calcFunc-integ (nth 1 expr) var low high)))
  1372.    ((and (eq (car-safe expr) '*)
  1373.      (not (math-expr-contains (nth 1 expr) var)))
  1374.     (math-mul (nth 1 expr) (calcFunc-integ (nth 2 expr) var low high)))
  1375.    ((and (eq (car-safe expr) '*)
  1376.      (not (math-expr-contains (nth 2 expr) var)))
  1377.     (math-mul (calcFunc-integ (nth 1 expr) var low high) (nth 2 expr)))
  1378.    ((and (eq (car-safe expr) '/)
  1379.      (not (math-expr-contains (nth 1 expr) var))
  1380.      (not (math-equal-int (nth 1 expr) 1)))
  1381.     (math-mul (nth 1 expr)
  1382.           (calcFunc-integ (math-div 1 (nth 2 expr)) var low high)))
  1383.    ((and (eq (car-safe expr) '/)
  1384.      (not (math-expr-contains (nth 2 expr) var)))
  1385.     (math-div (calcFunc-integ (nth 1 expr) var low high) (nth 2 expr)))
  1386.    ((and (eq (car-safe expr) '/)
  1387.      (eq (car-safe (nth 1 expr)) '*)
  1388.      (not (math-expr-contains (nth 1 (nth 1 expr)) var)))
  1389.     (math-mul (nth 1 (nth 1 expr))
  1390.           (calcFunc-integ (math-div (nth 2 (nth 1 expr)) (nth 2 expr))
  1391.                   var low high)))
  1392.    ((and (eq (car-safe expr) '/)
  1393.      (eq (car-safe (nth 1 expr)) '*)
  1394.      (not (math-expr-contains (nth 2 (nth 1 expr)) var)))
  1395.     (math-mul (nth 2 (nth 1 expr))
  1396.           (calcFunc-integ (math-div (nth 1 (nth 1 expr)) (nth 2 expr))
  1397.                   var low high)))
  1398.    ((and (eq (car-safe expr) '/)
  1399.      (eq (car-safe (nth 2 expr)) '*)
  1400.      (not (math-expr-contains (nth 1 (nth 2 expr)) var)))
  1401.     (math-div (calcFunc-integ (math-div (nth 1 expr) (nth 2 (nth 2 expr)))
  1402.                   var low high)
  1403.           (nth 1 (nth 2 expr))))
  1404.    ((and (eq (car-safe expr) '/)
  1405.      (eq (car-safe (nth 2 expr)) '*)
  1406.      (not (math-expr-contains (nth 2 (nth 2 expr)) var)))
  1407.     (math-div (calcFunc-integ (math-div (nth 1 expr) (nth 1 (nth 2 expr)))
  1408.                   var low high)
  1409.           (nth 2 (nth 2 expr))))
  1410.    ((eq (car-safe expr) 'vec)
  1411.     (cons 'vec (mapcar (function (lambda (x) (calcFunc-integ x var low high)))
  1412.                (cdr expr))))
  1413.    (t
  1414.     (let ((state (list calc-angle-mode
  1415.                ;;calc-symbolic-mode
  1416.                ;;calc-prefer-frac
  1417.                calc-internal-prec
  1418.                (calc-var-value 'var-IntegRules)
  1419.                (calc-var-value 'var-IntegSimpRules))))
  1420.       (or (equal state math-integral-cache-state)
  1421.       (setq math-integral-cache-state state
  1422.         math-integral-cache nil)))
  1423.     (let* ((math-max-integral-limit (or (and (boundp 'var-IntegLimit)
  1424.                          (natnump var-IntegLimit)
  1425.                          var-IntegLimit)
  1426.                     3))
  1427.        (math-integral-limit 1)
  1428.        (sexpr (math-expr-subst expr var math-integ-var))
  1429.        (trace-buffer (get-buffer "*Trace*"))
  1430.        (calc-language (if (eq calc-language 'big) nil calc-language))
  1431.        (math-any-substs t)
  1432.        (math-enable-subst nil)
  1433.        (math-prev-parts-v nil)
  1434.        (math-doing-parts nil)
  1435.        (math-good-parts nil)
  1436.        (res
  1437.         (if trace-buffer
  1438.         (let ((calcbuf (current-buffer))
  1439.               (calcwin (selected-window)))
  1440.           (unwind-protect
  1441.               (progn
  1442.             (if (get-buffer-window trace-buffer)
  1443.                 (select-window (get-buffer-window trace-buffer)))
  1444.             (set-buffer trace-buffer)
  1445.             (goto-char (point-max))
  1446.             (or (assq 'scroll-stop (buffer-local-variables))
  1447.                 (progn
  1448.                   (make-local-variable 'scroll-step)
  1449.                   (setq scroll-step 3)))
  1450.             (insert "\n\n\n")
  1451.             (set-buffer calcbuf)
  1452.             (math-try-integral sexpr))
  1453.             (select-window calcwin)
  1454.               (set-buffer calcbuf)))
  1455.           (math-try-integral sexpr))))
  1456.       (if res
  1457.       (progn
  1458.         (if (calc-has-rules 'var-IntegAfterRules)
  1459.         (setq res (math-rewrite res '(var IntegAfterRules
  1460.                           var-IntegAfterRules))))
  1461.         (math-simplify
  1462.          (if (and low high)
  1463.          (math-sub (math-expr-subst res math-integ-var high)
  1464.                (math-expr-subst res math-integ-var low))
  1465.            (setq res (math-fix-const-terms res math-integ-vars))
  1466.            (if low
  1467.            (math-expr-subst res math-integ-var low)
  1468.          (math-expr-subst res math-integ-var var)))))
  1469.     (append (list 'calcFunc-integ expr var)
  1470.         (and low (list low))
  1471.         (and high (list high)))))))
  1472. )
  1473.  
  1474.  
  1475. (math-defintegral calcFunc-inv
  1476.   (math-integral (math-div 1 u)))
  1477.  
  1478. (math-defintegral calcFunc-conj
  1479.   (let ((int (math-integral u)))
  1480.     (and int
  1481.      (list 'calcFunc-conj int))))
  1482.  
  1483. (math-defintegral calcFunc-deg
  1484.   (let ((int (math-integral u)))
  1485.     (and int
  1486.      (list 'calcFunc-deg int))))
  1487.  
  1488. (math-defintegral calcFunc-rad
  1489.   (let ((int (math-integral u)))
  1490.     (and int
  1491.      (list 'calcFunc-rad int))))
  1492.  
  1493. (math-defintegral calcFunc-re
  1494.   (let ((int (math-integral u)))
  1495.     (and int
  1496.      (list 'calcFunc-re int))))
  1497.  
  1498. (math-defintegral calcFunc-im
  1499.   (let ((int (math-integral u)))
  1500.     (and int
  1501.      (list 'calcFunc-im int))))
  1502.  
  1503. (math-defintegral calcFunc-sqrt
  1504.   (and (equal u math-integ-var)
  1505.        (math-mul '(frac 2 3)
  1506.          (list 'calcFunc-sqrt (math-pow u 3)))))
  1507.  
  1508. (math-defintegral calcFunc-exp
  1509.   (or (and (equal u math-integ-var)
  1510.        (list 'calcFunc-exp u))
  1511.       (let ((p (math-is-polynomial u math-integ-var 2)))
  1512.     (and (nth 2 p)
  1513.          (let ((sqa (math-sqrt (math-neg (nth 2 p)))))
  1514.            (math-div
  1515.         (math-mul
  1516.          (math-mul (math-div (list 'calcFunc-sqrt '(var pi var-pi))
  1517.                      sqa)
  1518.                (math-normalize
  1519.                 (list 'calcFunc-exp
  1520.                   (math-div (math-sub (math-mul (car p)
  1521.                                 (nth 2 p))
  1522.                               (math-div
  1523.                                (math-sqr (nth 1 p))
  1524.                                4))
  1525.                         (nth 2 p)))))
  1526.          (list 'calcFunc-erf
  1527.                (math-sub (math-mul sqa math-integ-var)
  1528.                  (math-div (nth 1 p) (math-mul 2 sqa)))))
  1529.         2))))))
  1530.  
  1531. (math-defintegral calcFunc-ln
  1532.   (or (and (equal u math-integ-var)
  1533.        (math-sub (math-mul u (list 'calcFunc-ln u)) u))
  1534.       (and (eq (car u) '*)
  1535.        (math-integral (math-add (list 'calcFunc-ln (nth 1 u))
  1536.                     (list 'calcFunc-ln (nth 2 u)))))
  1537.       (and (eq (car u) '/)
  1538.        (math-integral (math-sub (list 'calcFunc-ln (nth 1 u))
  1539.                     (list 'calcFunc-ln (nth 2 u)))))
  1540.       (and (eq (car u) '^)
  1541.        (math-integral (math-mul (nth 2 u)
  1542.                     (list 'calcFunc-ln (nth 1 u)))))))
  1543.  
  1544. (math-defintegral calcFunc-log10
  1545.   (and (equal u math-integ-var)
  1546.        (math-sub (math-mul u (list 'calcFunc-ln u))
  1547.          (math-div u (list 'calcFunc-ln 10)))))
  1548.  
  1549. (math-defintegral-2 calcFunc-log
  1550.   (math-integral (math-div (list 'calcFunc-ln u)
  1551.                (list 'calcFunc-ln v))))
  1552.  
  1553. (math-defintegral calcFunc-sin
  1554.   (or (and (equal u math-integ-var)
  1555.        (math-neg (math-from-radians-2 (list 'calcFunc-cos u))))
  1556.       (and (nth 2 (math-is-polynomial u math-integ-var 2))
  1557.        (math-integral (math-to-exponentials (list 'calcFunc-sin u))))))
  1558.  
  1559. (math-defintegral calcFunc-cos
  1560.   (or (and (equal u math-integ-var)
  1561.        (math-from-radians-2 (list 'calcFunc-sin u)))
  1562.       (and (nth 2 (math-is-polynomial u math-integ-var 2))
  1563.        (math-integral (math-to-exponentials (list 'calcFunc-cos u))))))
  1564.  
  1565. (math-defintegral calcFunc-tan
  1566.   (and (equal u math-integ-var)
  1567.        (math-neg (math-from-radians-2
  1568.           (list 'calcFunc-ln (list 'calcFunc-cos u))))))
  1569.  
  1570. (math-defintegral calcFunc-arcsin
  1571.   (and (equal u math-integ-var)
  1572.        (math-add (math-mul u (list 'calcFunc-arcsin u))
  1573.          (math-from-radians-2
  1574.           (list 'calcFunc-sqrt (math-sub 1 (math-sqr u)))))))
  1575.  
  1576. (math-defintegral calcFunc-arccos
  1577.   (and (equal u math-integ-var)
  1578.        (math-sub (math-mul u (list 'calcFunc-arccos u))
  1579.          (math-from-radians-2
  1580.           (list 'calcFunc-sqrt (math-sub 1 (math-sqr u)))))))
  1581.  
  1582. (math-defintegral calcFunc-arctan
  1583.   (and (equal u math-integ-var)
  1584.        (math-sub (math-mul u (list 'calcFunc-arctan u))
  1585.          (math-from-radians-2
  1586.           (math-div (list 'calcFunc-ln (math-add 1 (math-sqr u)))
  1587.                 2)))))
  1588.  
  1589. (math-defintegral calcFunc-sinh
  1590.   (and (equal u math-integ-var)
  1591.        (list 'calcFunc-cosh u)))
  1592.  
  1593. (math-defintegral calcFunc-cosh
  1594.   (and (equal u math-integ-var)
  1595.        (list 'calcFunc-sinh u)))
  1596.  
  1597. (math-defintegral calcFunc-tanh
  1598.   (and (equal u math-integ-var)
  1599.        (list 'calcFunc-ln (list 'calcFunc-cosh u))))
  1600.  
  1601. (math-defintegral calcFunc-arcsinh
  1602.   (and (equal u math-integ-var)
  1603.        (math-sub (math-mul u (list 'calcFunc-arcsinh u))
  1604.          (list 'calcFunc-sqrt (math-add (math-sqr u) 1)))))
  1605.  
  1606. (math-defintegral calcFunc-arccosh
  1607.   (and (equal u math-integ-var)
  1608.        (math-sub (math-mul u (list 'calcFunc-arccosh u))
  1609.          (list 'calcFunc-sqrt (math-sub 1 (math-sqr u))))))
  1610.  
  1611. (math-defintegral calcFunc-arctanh
  1612.   (and (equal u math-integ-var)
  1613.        (math-sub (math-mul u (list 'calcFunc-arctan u))
  1614.          (math-div (list 'calcFunc-ln
  1615.                  (math-add 1 (math-sqr u)))
  1616.                2))))
  1617.  
  1618. ;;; (Ax + B) / (ax^2 + bx + c)^n forms.
  1619. (math-defintegral-2 /
  1620.   (math-integral-rational-funcs u v))
  1621.  
  1622. (defun math-integral-rational-funcs (u v)
  1623.   (let ((pu (math-is-polynomial u math-integ-var 1))
  1624.     (vpow 1) pv)
  1625.     (and pu
  1626.      (catch 'int-rat
  1627.        (if (and (eq (car-safe v) '^) (natnump (nth 2 v)))
  1628.            (setq vpow (nth 2 v)
  1629.              v (nth 1 v)))
  1630.        (and (setq pv (math-is-polynomial v math-integ-var 2))
  1631.         (let ((int (math-mul-thru
  1632.                 (car pu)
  1633.                 (math-integral-q02 (car pv) (nth 1 pv)
  1634.                            (nth 2 pv) v vpow))))
  1635.           (if (cdr pu)
  1636.               (setq int (math-add int
  1637.                       (math-mul-thru
  1638.                        (nth 1 pu)
  1639.                        (math-integral-q12
  1640.                         (car pv) (nth 1 pv)
  1641.                         (nth 2 pv) v vpow)))))
  1642.           int))))))
  1643.  
  1644. (defun math-integral-q12 (a b c v vpow)
  1645.   (let (q)
  1646.     (cond ((not c)
  1647.        (cond ((= vpow 1)
  1648.           (math-sub (math-div math-integ-var b)
  1649.                 (math-mul (math-div a (math-sqr b))
  1650.                       (list 'calcFunc-ln v))))
  1651.          ((= vpow 2)
  1652.           (math-div (math-add (list 'calcFunc-ln v)
  1653.                       (math-div a v))
  1654.                 (math-sqr b)))
  1655.          (t
  1656.           (let ((nm1 (math-sub vpow 1))
  1657.             (nm2 (math-sub vpow 2)))
  1658.             (math-div (math-sub
  1659.                    (math-div a (math-mul nm1 (math-pow v nm1)))
  1660.                    (math-div 1 (math-mul nm2 (math-pow v nm2))))
  1661.                   (math-sqr b))))))
  1662.       ((math-zerop
  1663.         (setq q (math-sub (math-mul 4 (math-mul a c)) (math-sqr b))))
  1664.        (let ((part (math-div b (math-mul 2 c))))
  1665.          (math-mul-thru (math-pow c vpow)
  1666.                 (math-integral-q12 part 1 nil
  1667.                            (math-add math-integ-var part)
  1668.                            (* vpow 2)))))
  1669.       ((= vpow 1)
  1670.        (and (math-ratp q) (math-negp q)
  1671.         (let ((calc-symbolic-mode t))
  1672.           (math-ratp (math-sqrt (math-neg q))))
  1673.         (throw 'int-rat nil))  ; should have used calcFunc-apart first
  1674.        (math-sub (math-div (list 'calcFunc-ln v) (math-mul 2 c))
  1675.              (math-mul-thru (math-div b (math-mul 2 c))
  1676.                     (math-integral-q02 a b c v 1))))
  1677.       (t
  1678.        (let ((n (1- vpow)))
  1679.          (math-sub (math-neg (math-div
  1680.                   (math-add (math-mul b math-integ-var)
  1681.                         (math-mul 2 a))
  1682.                   (math-mul n (math-mul q (math-pow v n)))))
  1683.                (math-mul-thru (math-div (math-mul b (1- (* 2 n)))
  1684.                         (math-mul n q))
  1685.                       (math-integral-q02 a b c v n)))))))
  1686. )
  1687.  
  1688. (defun math-integral-q02 (a b c v vpow)
  1689.   (let (q rq part)
  1690.     (cond ((not c)
  1691.        (cond ((= vpow 1)
  1692.           (math-div (list 'calcFunc-ln v) b))
  1693.          (t
  1694.           (math-div (math-pow v (- 1 vpow))
  1695.                 (math-mul (- 1 vpow) b)))))
  1696.       ((math-zerop
  1697.         (setq q (math-sub (math-mul 4 (math-mul a c)) (math-sqr b))))
  1698.        (let ((part (math-div b (math-mul 2 c))))
  1699.          (math-mul-thru (math-pow c vpow)
  1700.                 (math-integral-q02 part 1 nil
  1701.                            (math-add math-integ-var part)
  1702.                            (* vpow 2)))))
  1703.       ((progn
  1704.          (setq part (math-add (math-mul 2 (math-mul c math-integ-var)) b))
  1705.          (> vpow 1))
  1706.        (let ((n (1- vpow)))
  1707.          (math-add (math-div part (math-mul n (math-mul q (math-pow v n))))
  1708.                (math-mul-thru (math-div (math-mul (- (* 4 n) 2) c)
  1709.                         (math-mul n q))
  1710.                       (math-integral-q02 a b c v n)))))
  1711.       ((math-guess-if-neg q)
  1712.        (setq rq (list 'calcFunc-sqrt (math-neg q)))
  1713.        ;;(math-div-thru (list 'calcFunc-ln
  1714.        ;;            (math-div (math-sub part rq)
  1715.        ;;                  (math-add part rq)))
  1716.        ;;          rq)
  1717.        (math-div (math-mul -2 (list 'calcFunc-arctanh
  1718.                     (math-div part rq)))
  1719.              rq))
  1720.       (t
  1721.        (setq rq (list 'calcFunc-sqrt q))
  1722.        (math-div (math-mul 2 (math-to-radians-2
  1723.                   (list 'calcFunc-arctan
  1724.                     (math-div part rq))))
  1725.              rq))))
  1726. )
  1727.  
  1728.  
  1729. (math-defintegral calcFunc-erf
  1730.   (and (equal u math-integ-var)
  1731.        (math-add (math-mul u (list 'calcFunc-erf u))
  1732.          (math-div 1 (math-mul (list 'calcFunc-exp (math-sqr u))
  1733.                        (list 'calcFunc-sqrt
  1734.                          '(var pi var-pi)))))))
  1735.  
  1736. (math-defintegral calcFunc-erfc
  1737.   (and (equal u math-integ-var)
  1738.        (math-sub (math-mul u (list 'calcFunc-erfc u))
  1739.          (math-div 1 (math-mul (list 'calcFunc-exp (math-sqr u))
  1740.                        (list 'calcFunc-sqrt
  1741.                          '(var pi var-pi)))))))
  1742.  
  1743.  
  1744.  
  1745.  
  1746. (defun calcFunc-table (expr var &optional low high step)
  1747.   (or low (setq low '(neg (var inf var-inf)) high '(var inf var-inf)))
  1748.   (or high (setq high low low 1))
  1749.   (and (or (math-infinitep low) (math-infinitep high))
  1750.        (not step)
  1751.        (math-scan-for-limits expr))
  1752.   (and step (math-zerop step) (math-reject-arg step 'nonzerop))
  1753.   (let ((known (+ (if (Math-objectp low) 1 0)
  1754.           (if (Math-objectp high) 1 0)
  1755.           (if (or (null step) (Math-objectp step)) 1 0)))
  1756.     (count '(var inf var-inf))
  1757.     vec)
  1758.     (or (= known 2)   ; handy optimization
  1759.     (equal high '(var inf var-inf))
  1760.     (progn
  1761.       (setq count (math-div (math-sub high low) (or step 1)))
  1762.       (or (Math-objectp count)
  1763.           (setq count (math-simplify count)))
  1764.       (if (Math-messy-integerp count)
  1765.           (setq count (math-trunc count)))))
  1766.     (if (Math-negp count)
  1767.     (setq count -1))
  1768.     (if (integerp count)
  1769.     (let ((var-DUMMY nil)
  1770.           (vec math-tabulate-initial)
  1771.           (math-working-step-2 (1+ count))
  1772.           (math-working-step 0))
  1773.       (setq expr (math-evaluate-expr
  1774.               (math-expr-subst expr var '(var DUMMY var-DUMMY))))
  1775.       (while (>= count 0)
  1776.         (setq math-working-step (1+ math-working-step)
  1777.           var-DUMMY low
  1778.           vec (cond ((eq math-tabulate-function 'calcFunc-sum)
  1779.                  (math-add vec (math-evaluate-expr expr)))
  1780.                 ((eq math-tabulate-function 'calcFunc-prod)
  1781.                  (math-mul vec (math-evaluate-expr expr)))
  1782.                 (t
  1783.                  (cons (math-evaluate-expr expr) vec)))
  1784.           low (math-add low (or step 1))
  1785.           count (1- count)))
  1786.       (if math-tabulate-function
  1787.           vec
  1788.         (cons 'vec (nreverse vec))))
  1789.       (if (Math-integerp count)
  1790.       (calc-record-why 'fixnump high)
  1791.     (if (Math-num-integerp low)
  1792.         (if (Math-num-integerp high)
  1793.         (calc-record-why 'integerp step)
  1794.           (calc-record-why 'integerp high))
  1795.       (calc-record-why 'integerp low)))
  1796.       (append (list (or math-tabulate-function 'calcFunc-table)
  1797.             expr var)
  1798.           (and (not (and (equal low '(neg (var inf var-inf)))
  1799.                  (equal high '(var inf var-inf))))
  1800.            (list low high))
  1801.           (and step (list step)))))
  1802. )
  1803.  
  1804. (setq math-tabulate-initial nil)
  1805. (setq math-tabulate-function nil)
  1806.  
  1807. (defun math-scan-for-limits (x)
  1808.   (cond ((Math-primp x))
  1809.     ((and (eq (car x) 'calcFunc-subscr)
  1810.           (Math-vectorp (nth 1 x))
  1811.           (math-expr-contains (nth 2 x) var))
  1812.      (let* ((calc-next-why nil)
  1813.         (low-val (math-solve-for (nth 2 x) 1 var nil))
  1814.         (high-val (math-solve-for (nth 2 x) (1- (length (nth 1 x)))
  1815.                       var nil))
  1816.         temp)
  1817.        (and low-val (math-realp low-val)
  1818.         high-val (math-realp high-val))
  1819.        (and (Math-lessp high-val low-val)
  1820.         (setq temp low-val low-val high-val high-val temp))
  1821.        (setq low (math-max low (math-ceiling low-val))
  1822.          high (math-min high (math-floor high-val)))))
  1823.     (t
  1824.      (while (setq x (cdr x))
  1825.        (math-scan-for-limits (car x)))))
  1826. )
  1827.  
  1828.  
  1829. (defun calcFunc-sum (expr var &optional low high step)
  1830.   (if math-disable-sums (math-reject-arg))
  1831.   (let* ((res (let* ((calc-internal-prec (+ calc-internal-prec 2)))
  1832.         (math-sum-rec expr var low high step)))
  1833.      (math-disable-sums t))
  1834.     (math-normalize res))
  1835. )
  1836. (setq math-disable-sums nil)
  1837.  
  1838. (defun math-sum-rec (expr var &optional low high step)
  1839.   (or low (setq low '(neg (var inf var-inf)) high '(var inf var-inf)))
  1840.   (and low (not high) (setq high low low 1))
  1841.   (let (t1 t2 val)
  1842.     (setq val
  1843.       (cond
  1844.        ((not (math-expr-contains expr var))
  1845.         (math-mul expr (math-add (math-div (math-sub high low) (or step 1))
  1846.                      1)))
  1847.        ((and step (not (math-equal-int step 1)))
  1848.         (if (math-negp step)
  1849.         (math-sum-rec expr var high low (math-neg step))
  1850.           (let ((lo (math-simplify (math-div low step))))
  1851.         (if (math-known-num-integerp lo)
  1852.             (math-sum-rec (math-normalize
  1853.                    (math-expr-subst expr var
  1854.                             (math-mul step var)))
  1855.                   var lo (math-simplify (math-div high step)))
  1856.           (math-sum-rec (math-normalize
  1857.                  (math-expr-subst expr var
  1858.                           (math-add (math-mul step var)
  1859.                                 low)))
  1860.                 var 0
  1861.                 (math-simplify (math-div (math-sub high low)
  1862.                              step)))))))
  1863.        ((memq (setq t1 (math-compare low high)) '(0 1))
  1864.         (if (eq t1 0)
  1865.         (math-expr-subst expr var low)
  1866.           0))
  1867.        ((setq t1 (math-is-polynomial expr var 20))
  1868.         (let ((poly nil)
  1869.           (n 0))
  1870.           (while t1
  1871.         (setq poly (math-poly-mix poly 1
  1872.                       (math-sum-integer-power n) (car t1))
  1873.               n (1+ n)
  1874.               t1 (cdr t1)))
  1875.           (setq n (math-build-polynomial-expr poly high))
  1876.           (if (memq low '(0 1))
  1877.           n
  1878.         (math-sub n (math-build-polynomial-expr poly
  1879.                             (math-sub low 1))))))
  1880.        ((and (memq (car expr) '(+ -))
  1881.          (setq t1 (math-sum-rec (nth 1 expr) var low high)
  1882.                t2 (math-sum-rec (nth 2 expr) var low high))
  1883.          (not (and (math-expr-calls t1 '(calcFunc-sum))
  1884.                (math-expr-calls t2 '(calcFunc-sum)))))
  1885.         (list (car expr) t1 t2))
  1886.        ((and (eq (car expr) '*)
  1887.          (setq t1 (math-sum-const-factors expr var)))
  1888.         (math-mul (car t1) (math-sum-rec (cdr t1) var low high)))
  1889.        ((and (eq (car expr) '*) (memq (car-safe (nth 1 expr)) '(+ -)))
  1890.         (math-sum-rec (math-add-or-sub (math-mul (nth 1 (nth 1 expr))
  1891.                              (nth 2 expr))
  1892.                        (math-mul (nth 2 (nth 1 expr))
  1893.                              (nth 2 expr))
  1894.                        nil (eq (car (nth 1 expr)) '-))
  1895.               var low high))
  1896.        ((and (eq (car expr) '*) (memq (car-safe (nth 2 expr)) '(+ -)))
  1897.         (math-sum-rec (math-add-or-sub (math-mul (nth 1 expr)
  1898.                              (nth 1 (nth 2 expr)))
  1899.                        (math-mul (nth 1 expr)
  1900.                              (nth 2 (nth 2 expr)))
  1901.                        nil (eq (car (nth 2 expr)) '-))
  1902.               var low high))
  1903.        ((and (eq (car expr) '/)
  1904.          (not (math-primp (nth 1 expr)))
  1905.          (setq t1 (math-sum-const-factors (nth 1 expr) var)))
  1906.         (math-mul (car t1)
  1907.               (math-sum-rec (math-div (cdr t1) (nth 2 expr))
  1908.                     var low high)))
  1909.        ((and (eq (car expr) '/)
  1910.          (setq t1 (math-sum-const-factors (nth 2 expr) var)))
  1911.         (math-div (math-sum-rec (math-div (nth 1 expr) (cdr t1))
  1912.                     var low high)
  1913.               (car t1)))
  1914.        ((eq (car expr) 'neg)
  1915.         (math-neg (math-sum-rec (nth 1 expr) var low high)))
  1916.        ((and (eq (car expr) '^)
  1917.          (not (math-expr-contains (nth 1 expr) var))
  1918.          (setq t1 (math-is-polynomial (nth 2 expr) var 1)))
  1919.         (let ((x (math-pow (nth 1 expr) (nth 1 t1))))
  1920.           (math-div (math-mul (math-sub (math-pow x (math-add 1 high))
  1921.                         (math-pow x low))
  1922.                   (math-pow (nth 1 expr) (car t1)))
  1923.             (math-sub x 1))))
  1924.        ((and (setq t1 (math-to-exponentials expr))
  1925.          (setq t1 (math-sum-rec t1 var low high))
  1926.          (not (math-expr-calls t1 '(calcFunc-sum))))
  1927.         (math-to-exps t1))
  1928.        ((memq (car expr) '(calcFunc-ln calcFunc-log10))
  1929.         (list (car expr) (calcFunc-prod (nth 1 expr) var low high)))
  1930.        ((and (eq (car expr) 'calcFunc-log)
  1931.          (= (length expr) 3)
  1932.          (not (math-expr-contains (nth 2 expr) var)))
  1933.         (list 'calcFunc-log
  1934.           (calcFunc-prod (nth 1 expr) var low high)
  1935.           (nth 2 expr)))))
  1936.     (if (equal val '(var nan var-nan)) (setq val nil))
  1937.     (or val
  1938.     (let* ((math-tabulate-initial 0)
  1939.            (math-tabulate-function 'calcFunc-sum))
  1940.       (calcFunc-table expr var low high))))
  1941. )
  1942.  
  1943. (defun calcFunc-asum (expr var low &optional high step no-mul-flag)
  1944.   (or high (setq high low low 1))
  1945.   (if (and step (not (math-equal-int step 1)))
  1946.       (if (math-negp step)
  1947.       (math-mul (math-pow -1 low)
  1948.             (calcFunc-asum expr var high low (math-neg step) t))
  1949.     (let ((lo (math-simplify (math-div low step))))
  1950.       (if (math-num-integerp lo)
  1951.           (calcFunc-asum (math-normalize
  1952.                   (math-expr-subst expr var
  1953.                            (math-mul step var)))
  1954.                  var lo (math-simplify (math-div high step)))
  1955.         (calcFunc-asum (math-normalize
  1956.                 (math-expr-subst expr var
  1957.                          (math-add (math-mul step var)
  1958.                                low)))
  1959.                var 0
  1960.                (math-simplify (math-div (math-sub high low)
  1961.                             step))))))
  1962.     (math-mul (if no-mul-flag 1 (math-pow -1 low))
  1963.           (calcFunc-sum (math-mul (math-pow -1 var) expr) var low high)))
  1964. )
  1965.  
  1966. (defun math-sum-const-factors (expr var)
  1967.   (let ((const nil)
  1968.     (not-const nil)
  1969.     (p expr))
  1970.     (while (eq (car-safe p) '*)
  1971.       (if (math-expr-contains (nth 1 p) var)
  1972.       (setq not-const (cons (nth 1 p) not-const))
  1973.     (setq const (cons (nth 1 p) const)))
  1974.       (setq p (nth 2 p)))
  1975.     (if (math-expr-contains p var)
  1976.     (setq not-const (cons p not-const))
  1977.       (setq const (cons p const)))
  1978.     (and const
  1979.      (cons (let ((temp (car const)))
  1980.          (while (setq const (cdr const))
  1981.            (setq temp (list '* (car const) temp)))
  1982.          temp)
  1983.            (let ((temp (or (car not-const) 1)))
  1984.          (while (setq not-const (cdr not-const))
  1985.            (setq temp (list '* (car not-const) temp)))
  1986.          temp))))
  1987. )
  1988.  
  1989. ;; Following is from CRC Math Tables, 27th ed, pp. 52-53.
  1990. (defun math-sum-integer-power (pow)
  1991.   (let ((calc-prefer-frac t)
  1992.     (n (length math-sum-int-pow-cache)))
  1993.     (while (<= n pow)
  1994.       (let* ((new (list 0 0))
  1995.          (lin new)
  1996.          (pp (cdr (nth (1- n) math-sum-int-pow-cache)))
  1997.          (p 2)
  1998.          (sum 0)
  1999.          q)
  2000.     (while pp
  2001.       (setq q (math-div (car pp) p)
  2002.         new (cons (math-mul q n) new)
  2003.         sum (math-add sum q)
  2004.         p (1+ p)
  2005.         pp (cdr pp)))
  2006.     (setcar lin (math-sub 1 (math-mul n sum)))
  2007.     (setq math-sum-int-pow-cache
  2008.           (nconc math-sum-int-pow-cache (list (nreverse new)))
  2009.           n (1+ n))))
  2010.     (nth pow math-sum-int-pow-cache))
  2011. )
  2012. (setq math-sum-int-pow-cache (list '(0 1)))
  2013.  
  2014. (defun math-to-exponentials (expr)
  2015.   (and (consp expr)
  2016.        (= (length expr) 2)
  2017.        (let ((x (nth 1 expr))
  2018.          (pi (if calc-symbolic-mode '(var pi var-pi) (math-pi)))
  2019.          (i (if calc-symbolic-mode '(var i var-i) '(cplx 0 1))))
  2020.      (cond ((eq (car expr) 'calcFunc-exp)
  2021.         (list '^ '(var e var-e) x))
  2022.            ((eq (car expr) 'calcFunc-sin)
  2023.         (or (eq calc-angle-mode 'rad)
  2024.             (setq x (list '/ (list '* x pi) 180)))
  2025.         (list '/ (list '-
  2026.                    (list '^ '(var e var-e) (list '* x i))
  2027.                    (list '^ '(var e var-e)
  2028.                      (list 'neg (list '* x i))))
  2029.               (list '* 2 i)))
  2030.            ((eq (car expr) 'calcFunc-cos)
  2031.         (or (eq calc-angle-mode 'rad)
  2032.             (setq x (list '/ (list '* x pi) 180)))
  2033.         (list '/ (list '+
  2034.                    (list '^ '(var e var-e)
  2035.                      (list '* x i))
  2036.                    (list '^ '(var e var-e)
  2037.                      (list 'neg (list '* x i))))
  2038.               2))
  2039.            ((eq (car expr) 'calcFunc-sinh)
  2040.         (list '/ (list '-
  2041.                    (list '^ '(var e var-e) x)
  2042.                    (list '^ '(var e var-e) (list 'neg x)))
  2043.               2))
  2044.            ((eq (car expr) 'calcFunc-cosh)
  2045.         (list '/ (list '+
  2046.                    (list '^ '(var e var-e) x)
  2047.                    (list '^ '(var e var-e) (list 'neg x)))
  2048.               2))
  2049.            (t nil))))
  2050. )
  2051.  
  2052. (defun math-to-exps (expr)
  2053.   (cond (calc-symbolic-mode expr)
  2054.     ((Math-primp expr)
  2055.      (if (equal expr '(var e var-e)) (math-e) expr))
  2056.     ((and (eq (car expr) '^)
  2057.           (equal (nth 1 expr) '(var e var-e)))
  2058.      (list 'calcFunc-exp (nth 2 expr)))
  2059.     (t
  2060.      (cons (car expr) (mapcar 'math-to-exps (cdr expr)))))
  2061. )
  2062.  
  2063.  
  2064. (defun calcFunc-prod (expr var &optional low high step)
  2065.   (if math-disable-prods (math-reject-arg))
  2066.   (let* ((res (let* ((calc-internal-prec (+ calc-internal-prec 2)))
  2067.         (math-prod-rec expr var low high step)))
  2068.      (math-disable-prods t))
  2069.     (math-normalize res))
  2070. )
  2071. (setq math-disable-prods nil)
  2072.  
  2073. (defun math-prod-rec (expr var &optional low high step)
  2074.   (or low (setq low '(neg (var inf var-inf)) high '(var inf var-inf)))
  2075.   (and low (not high) (setq high '(var inf var-inf)))
  2076.   (let (t1 t2 t3 val)
  2077.     (setq val
  2078.       (cond
  2079.        ((not (math-expr-contains expr var))
  2080.         (math-pow expr (math-add (math-div (math-sub high low) (or step 1))
  2081.                      1)))
  2082.        ((and step (not (math-equal-int step 1)))
  2083.         (if (math-negp step)
  2084.         (math-prod-rec expr var high low (math-neg step))
  2085.           (let ((lo (math-simplify (math-div low step))))
  2086.         (if (math-known-num-integerp lo)
  2087.             (math-prod-rec (math-normalize
  2088.                     (math-expr-subst expr var
  2089.                              (math-mul step var)))
  2090.                    var lo (math-simplify (math-div high step)))
  2091.           (math-prod-rec (math-normalize
  2092.                   (math-expr-subst expr var
  2093.                            (math-add (math-mul step
  2094.                                        var)
  2095.                                  low)))
  2096.                  var 0
  2097.                  (math-simplify (math-div (math-sub high low)
  2098.                               step)))))))
  2099.        ((and (memq (car expr) '(* /))
  2100.          (setq t1 (math-prod-rec (nth 1 expr) var low high)
  2101.                t2 (math-prod-rec (nth 2 expr) var low high))
  2102.          (not (and (math-expr-calls t1 '(calcFunc-prod))
  2103.                (math-expr-calls t2 '(calcFunc-prod)))))
  2104.         (list (car expr) t1 t2))
  2105.        ((and (eq (car expr) '^)
  2106.          (not (math-expr-contains (nth 2 expr) var)))
  2107.         (math-pow (math-prod-rec (nth 1 expr) var low high)
  2108.               (nth 2 expr)))
  2109.        ((and (eq (car expr) '^)
  2110.          (not (math-expr-contains (nth 1 expr) var)))
  2111.         (math-pow (nth 1 expr)
  2112.               (calcFunc-sum (nth 2 expr) var low high)))
  2113.        ((eq (car expr) 'sqrt)
  2114.         (math-normalize (list 'calcFunc-sqrt
  2115.                   (list 'calcFunc-prod (nth 1 expr)
  2116.                     var low high))))
  2117.        ((eq (car expr) 'neg)
  2118.         (math-mul (math-pow -1 (math-add (math-sub high low) 1))
  2119.               (math-prod-rec (nth 1 expr) var low high)))
  2120.        ((eq (car expr) 'calcFunc-exp)
  2121.         (list 'calcFunc-exp (calcFunc-sum (nth 1 expr) var low high)))
  2122.        ((and (setq t1 (math-is-polynomial expr var 1))
  2123.          (setq t2
  2124.                (cond
  2125.             ((or (and (math-equal-int (nth 1 t1) 1)
  2126.                   (setq low (math-simplify
  2127.                          (math-add low (car t1)))
  2128.                     high (math-simplify
  2129.                           (math-add high (car t1)))))
  2130.                  (and (math-equal-int (nth 1 t1) -1)
  2131.                   (setq t2 low
  2132.                     low (math-simplify
  2133.                          (math-sub (car t1) high))
  2134.                     high (math-simplify
  2135.                           (math-sub (car t1) t2)))))
  2136.              (if (or (math-zerop low) (math-zerop high))
  2137.                  0
  2138.                (if (and (or (math-negp low) (math-negp high))
  2139.                     (or (math-num-integerp low)
  2140.                     (math-num-integerp high)))
  2141.                    (if (math-posp high)
  2142.                    0
  2143.                  (math-mul (math-pow -1
  2144.                              (math-add
  2145.                               (math-add low high) 1))
  2146.                        (list '/
  2147.                          (list 'calcFunc-fact
  2148.                                (math-neg low))
  2149.                          (list 'calcFunc-fact
  2150.                                (math-sub -1 high)))))
  2151.                  (list '/
  2152.                    (list 'calcFunc-fact high)
  2153.                    (list 'calcFunc-fact (math-sub low 1))))))
  2154.             ((and (or (and (math-equal-int (nth 1 t1) 2)
  2155.                        (setq t2 (math-simplify
  2156.                          (math-add (math-mul low 2)
  2157.                                (car t1)))
  2158.                          t3 (math-simplify
  2159.                          (math-add (math-mul high 2)
  2160.                                (car t1)))))
  2161.                   (and (math-equal-int (nth 1 t1) -2)
  2162.                        (setq t2 (math-simplify
  2163.                          (math-sub (car t1)
  2164.                                (math-mul high 2)))
  2165.                          t3 (math-simplify 
  2166.                          (math-sub (car t1)
  2167.                                (math-mul low
  2168.                                      2))))))
  2169.                   (or (math-integerp t2)
  2170.                   (and (math-messy-integerp t2)
  2171.                        (setq t2 (math-trunc t2)))
  2172.                   (math-integerp t3)
  2173.                   (and (math-messy-integerp t3)
  2174.                        (setq t3 (math-trunc t3)))))
  2175.              (if (or (math-zerop t2) (math-zerop t3))
  2176.                  0
  2177.                (if (or (math-evenp t2) (math-evenp t3))
  2178.                    (if (or (math-negp t2) (math-negp t3))
  2179.                    (if (math-posp high)
  2180.                        0
  2181.                      (list '/
  2182.                        (list 'calcFunc-dfact
  2183.                          (math-neg t2))
  2184.                        (list 'calcFunc-dfact
  2185.                          (math-sub -2 t3))))
  2186.                  (list '/
  2187.                        (list 'calcFunc-dfact t3)
  2188.                        (list 'calcFunc-dfact
  2189.                          (math-sub t2 2))))
  2190.                  (if (math-negp t3)
  2191.                  (list '*
  2192.                        (list '^ -1
  2193.                          (list '/ (list '- (list '- t2 t3)
  2194.                                 2)
  2195.                            2))
  2196.                        (list '/
  2197.                          (list 'calcFunc-dfact
  2198.                            (math-neg t2))
  2199.                          (list 'calcFunc-dfact
  2200.                            (math-sub -2 t3))))
  2201.                    (if (math-posp t2)
  2202.                    (list '/
  2203.                      (list 'calcFunc-dfact t3)
  2204.                      (list 'calcFunc-dfact
  2205.                            (math-sub t2 2)))
  2206.                  nil))))))))
  2207.         t2)))
  2208.     (if (equal val '(var nan var-nan)) (setq val nil))
  2209.     (or val
  2210.     (let* ((math-tabulate-initial 1)
  2211.            (math-tabulate-function 'calcFunc-prod))
  2212.       (calcFunc-table expr var low high))))
  2213. )
  2214.  
  2215.  
  2216.  
  2217.  
  2218. ;;; Attempt to reduce lhs = rhs to solve-var = rhs', where solve-var appears
  2219. ;;; in lhs but not in rhs or rhs'; return rhs'.
  2220. ;;; Uses global values: solve-*.
  2221. (defun math-try-solve-for (lhs rhs &optional sign no-poly)
  2222.   (let (t1 t2 t3)
  2223.     (cond ((equal lhs solve-var)
  2224.        (setq math-solve-sign sign)
  2225.        (if (eq solve-full 'all)
  2226.            (let ((vec (list 'vec (math-evaluate-expr rhs)))
  2227.              newvec var p)
  2228.          (while math-solve-ranges
  2229.            (setq p (car math-solve-ranges)
  2230.              var (car p)
  2231.              newvec (list 'vec))
  2232.            (while (setq p (cdr p))
  2233.              (setq newvec (nconc newvec
  2234.                      (cdr (math-expr-subst
  2235.                            vec var (car p))))))
  2236.            (setq vec newvec
  2237.              math-solve-ranges (cdr math-solve-ranges)))
  2238.          (math-normalize vec))
  2239.          rhs))
  2240.       ((Math-primp lhs)
  2241.        nil)
  2242.       ((and (eq (car lhs) '-)
  2243.         (eq (car-safe (nth 1 lhs)) (car-safe (nth 2 lhs)))
  2244.         (Math-zerop rhs)
  2245.         (= (length (nth 1 lhs)) 2)
  2246.         (= (length (nth 2 lhs)) 2)
  2247.         (setq t1 (get (car (nth 1 lhs)) 'math-inverse))
  2248.         (setq t2 (funcall t1 '(var SOLVEDUM SOLVEDUM)))
  2249.         (eq (math-expr-contains-count t2 '(var SOLVEDUM SOLVEDUM)) 1)
  2250.         (setq t3 (math-solve-above-dummy t2))
  2251.         (setq t1 (math-try-solve-for (math-sub (nth 1 (nth 1 lhs))
  2252.                                (math-expr-subst
  2253.                             t2 t3
  2254.                             (nth 1 (nth 2 lhs))))
  2255.                          0)))
  2256.        t1)
  2257.       ((eq (car lhs) 'neg)
  2258.        (math-try-solve-for (nth 1 lhs) (math-neg rhs)
  2259.                    (and sign (- sign))))
  2260.       ((and (not (eq solve-full 't)) (math-try-solve-prod)))
  2261.       ((and (not no-poly)
  2262.         (setq t2 (math-decompose-poly lhs solve-var 15 rhs)))
  2263.        (setq t1 (cdr (nth 1 t2))
  2264.          t1 (let ((math-solve-ranges math-solve-ranges))
  2265.               (cond ((= (length t1) 5)
  2266.                  (apply 'math-solve-quartic (car t2) t1))
  2267.                 ((= (length t1) 4)
  2268.                  (apply 'math-solve-cubic (car t2) t1))
  2269.                 ((= (length t1) 3)
  2270.                  (apply 'math-solve-quadratic (car t2) t1))
  2271.                 ((= (length t1) 2)
  2272.                  (apply 'math-solve-linear (car t2) sign t1))
  2273.                 (solve-full
  2274.                  (math-poly-all-roots (car t2) t1))
  2275.                 (calc-symbolic-mode nil)
  2276.                 (t
  2277.                  (math-try-solve-for
  2278.                   (car t2)
  2279.                   (math-poly-any-root (reverse t1) 0 t)
  2280.                   nil t)))))
  2281.        (if t1
  2282.            (if (eq (nth 2 t2) 1)
  2283.            t1
  2284.          (math-solve-prod t1 (math-try-solve-for (nth 2 t2) 0 nil t)))
  2285.          (calc-record-why "*Unable to find a symbolic solution")
  2286.          nil))
  2287.       ((and (math-solve-find-root-term lhs nil)
  2288.         (eq (math-expr-contains-count lhs t1) 1))   ; just in case
  2289.        (math-try-solve-for (math-simplify
  2290.                 (math-sub (if (or t3 (math-evenp t2))
  2291.                           (math-pow t1 t2)
  2292.                         (math-neg (math-pow t1 t2)))
  2293.                       (math-expand-power
  2294.                        (math-sub (math-normalize
  2295.                               (math-expr-subst
  2296.                                lhs t1 0))
  2297.                              rhs)
  2298.                        t2 solve-var)))
  2299.                    0))
  2300.       ((eq (car lhs) '+)
  2301.        (cond ((not (math-expr-contains (nth 1 lhs) solve-var))
  2302.           (math-try-solve-for (nth 2 lhs)
  2303.                       (math-sub rhs (nth 1 lhs))
  2304.                       sign))
  2305.          ((not (math-expr-contains (nth 2 lhs) solve-var))
  2306.           (math-try-solve-for (nth 1 lhs)
  2307.                       (math-sub rhs (nth 2 lhs))
  2308.                       sign))))
  2309.       ((eq (car lhs) 'calcFunc-eq)
  2310.        (math-try-solve-for (math-sub (nth 1 lhs) (nth 2 lhs))
  2311.                    rhs sign no-poly))
  2312.       ((eq (car lhs) '-)
  2313.        (cond ((or (and (eq (car-safe (nth 1 lhs)) 'calcFunc-sin)
  2314.                (eq (car-safe (nth 2 lhs)) 'calcFunc-cos))
  2315.               (and (eq (car-safe (nth 1 lhs)) 'calcFunc-cos)
  2316.                (eq (car-safe (nth 2 lhs)) 'calcFunc-sin)))
  2317.           (math-try-solve-for (math-sub (nth 1 lhs)
  2318.                         (list (car (nth 1 lhs))
  2319.                               (math-sub
  2320.                                (math-quarter-circle t)
  2321.                                (nth 1 (nth 2 lhs)))))
  2322.                       rhs))
  2323.          ((not (math-expr-contains (nth 1 lhs) solve-var))
  2324.           (math-try-solve-for (nth 2 lhs)
  2325.                       (math-sub (nth 1 lhs) rhs)
  2326.                       (and sign (- sign))))
  2327.          ((not (math-expr-contains (nth 2 lhs) solve-var))
  2328.           (math-try-solve-for (nth 1 lhs)
  2329.                       (math-add rhs (nth 2 lhs))
  2330.                       sign))))
  2331.       ((and (eq solve-full 't) (math-try-solve-prod)))
  2332.       ((and (eq (car lhs) '%)
  2333.         (not (math-expr-contains (nth 2 lhs) solve-var)))
  2334.        (math-try-solve-for (nth 1 lhs) (math-add rhs
  2335.                              (math-solve-get-int
  2336.                               (nth 2 lhs)))))
  2337.       ((eq (car lhs) 'calcFunc-log)
  2338.        (cond ((not (math-expr-contains (nth 2 lhs) solve-var))
  2339.           (math-try-solve-for (nth 1 lhs) (math-pow (nth 2 lhs) rhs)))
  2340.          ((not (math-expr-contains (nth 1 lhs) solve-var))
  2341.           (math-try-solve-for (nth 2 lhs) (math-pow
  2342.                            (nth 1 lhs)
  2343.                            (math-div 1 rhs))))))
  2344.       ((and (= (length lhs) 2)
  2345.         (symbolp (car lhs))
  2346.         (setq t1 (get (car lhs) 'math-inverse))
  2347.         (setq t2 (funcall t1 rhs)))
  2348.        (setq t1 (get (car lhs) 'math-inverse-sign))
  2349.        (math-try-solve-for (nth 1 lhs) (math-normalize t2)
  2350.                    (and sign t1
  2351.                     (if (integerp t1)
  2352.                     (* t1 sign)
  2353.                       (funcall t1 lhs sign)))))
  2354.       ((and (symbolp (car lhs))
  2355.         (setq t1 (get (car lhs) 'math-inverse-n))
  2356.         (setq t2 (funcall t1 lhs rhs)))
  2357.        t2)
  2358.       ((setq t1 (math-expand-formula lhs))
  2359.        (math-try-solve-for t1 rhs sign))
  2360.       (t
  2361.        (calc-record-why "*No inverse known" lhs)
  2362.        nil)))
  2363. )
  2364.  
  2365. (setq math-solve-ranges nil)
  2366.  
  2367. (defun math-try-solve-prod ()
  2368.   (cond ((eq (car lhs) '*)
  2369.      (cond ((not (math-expr-contains (nth 1 lhs) solve-var))
  2370.         (math-try-solve-for (nth 2 lhs)
  2371.                     (math-div rhs (nth 1 lhs))
  2372.                     (math-solve-sign sign (nth 1 lhs))))
  2373.            ((not (math-expr-contains (nth 2 lhs) solve-var))
  2374.         (math-try-solve-for (nth 1 lhs)
  2375.                     (math-div rhs (nth 2 lhs))
  2376.                     (math-solve-sign sign (nth 2 lhs))))
  2377.            ((Math-zerop rhs)
  2378.         (math-solve-prod (let ((math-solve-ranges math-solve-ranges))
  2379.                    (math-try-solve-for (nth 2 lhs) 0))
  2380.                  (math-try-solve-for (nth 1 lhs) 0)))))
  2381.     ((eq (car lhs) '/)
  2382.      (cond ((not (math-expr-contains (nth 1 lhs) solve-var))
  2383.         (math-try-solve-for (nth 2 lhs)
  2384.                     (math-div (nth 1 lhs) rhs)
  2385.                     (math-solve-sign sign (nth 1 lhs))))
  2386.            ((not (math-expr-contains (nth 2 lhs) solve-var))
  2387.         (math-try-solve-for (nth 1 lhs)
  2388.                     (math-mul rhs (nth 2 lhs))
  2389.                     (math-solve-sign sign (nth 2 lhs))))
  2390.            ((setq t1 (math-try-solve-for (math-sub (nth 1 lhs)
  2391.                                (math-mul (nth 2 lhs)
  2392.                                  rhs))
  2393.                          0))
  2394.         t1)))
  2395.     ((eq (car lhs) '^)
  2396.      (cond ((not (math-expr-contains (nth 1 lhs) solve-var))
  2397.         (math-try-solve-for
  2398.          (nth 2 lhs)
  2399.          (math-add (math-normalize
  2400.                 (list 'calcFunc-log rhs (nth 1 lhs)))
  2401.                (math-div
  2402.                 (math-mul 2
  2403.                       (math-mul '(var pi var-pi)
  2404.                         (math-solve-get-int
  2405.                          '(var i var-i))))
  2406.                 (math-normalize
  2407.                  (list 'calcFunc-ln (nth 1 lhs)))))))
  2408.            ((not (math-expr-contains (nth 2 lhs) solve-var))
  2409.         (cond ((and (integerp (nth 2 lhs))
  2410.                 (>= (nth 2 lhs) 2)
  2411.                 (setq t1 (math-integer-log2 (nth 2 lhs))))
  2412.                (setq t2 rhs)
  2413.                (if (and (eq solve-full t)
  2414.                 (math-known-realp (nth 1 lhs)))
  2415.                (progn
  2416.                  (while (>= (setq t1 (1- t1)) 0)
  2417.                    (setq t2 (list 'calcFunc-sqrt t2)))
  2418.                  (setq t2 (math-solve-get-sign t2)))
  2419.              (while (>= (setq t1 (1- t1)) 0)
  2420.                (setq t2 (math-solve-get-sign
  2421.                      (math-normalize
  2422.                       (list 'calcFunc-sqrt t2))))))
  2423.                (math-try-solve-for
  2424.             (nth 1 lhs)
  2425.             (math-normalize t2)))
  2426.               ((math-looks-negp (nth 2 lhs))
  2427.                (math-try-solve-for
  2428.             (list '^ (nth 1 lhs) (math-neg (nth 2 lhs)))
  2429.             (math-div 1 rhs)))
  2430.               ((and (eq solve-full t)
  2431.                 (Math-integerp (nth 2 lhs))
  2432.                 (math-known-realp (nth 1 lhs)))
  2433.                (setq t1 (math-normalize
  2434.                  (list 'calcFunc-nroot rhs (nth 2 lhs))))
  2435.                (if (math-evenp (nth 2 lhs))
  2436.                (setq t1 (math-solve-get-sign t1)))
  2437.                (math-try-solve-for
  2438.             (nth 1 lhs) t1
  2439.             (and sign
  2440.                  (math-oddp (nth 2 lhs))
  2441.                  (math-solve-sign sign (nth 2 lhs)))))
  2442.               (t (math-try-solve-for
  2443.               (nth 1 lhs)
  2444.               (math-mul
  2445.                (math-normalize
  2446.                 (list 'calcFunc-exp
  2447.                   (if (Math-realp (nth 2 lhs))
  2448.                       (math-div (math-mul
  2449.                          '(var pi var-pi)
  2450.                          (math-solve-get-int
  2451.                           '(var i var-i)
  2452.                           (and (integerp (nth 2 lhs))
  2453.                                (math-abs
  2454.                             (nth 2 lhs)))))
  2455.                         (math-div (nth 2 lhs) 2))
  2456.                     (math-div (math-mul
  2457.                            2
  2458.                            (math-mul
  2459.                         '(var pi var-pi)
  2460.                         (math-solve-get-int
  2461.                          '(var i var-i)
  2462.                          (and (integerp (nth 2 lhs))
  2463.                               (math-abs
  2464.                                (nth 2 lhs))))))
  2465.                           (nth 2 lhs)))))
  2466.                (math-normalize
  2467.                 (list 'calcFunc-nroot
  2468.                   rhs
  2469.                   (nth 2 lhs))))
  2470.               (and sign
  2471.                    (math-oddp (nth 2 lhs))
  2472.                    (math-solve-sign sign (nth 2 lhs)))))))))
  2473.     (t nil))
  2474. )
  2475.  
  2476. (defun math-solve-prod (lsoln rsoln)
  2477.   (cond ((null lsoln)
  2478.      rsoln)
  2479.     ((null rsoln)
  2480.      lsoln)
  2481.     ((eq solve-full 'all)
  2482.      (cons 'vec (append (cdr lsoln) (cdr rsoln))))
  2483.     (solve-full
  2484.      (list 'calcFunc-if
  2485.            (list 'calcFunc-gt (math-solve-get-sign 1) 0)
  2486.            lsoln
  2487.            rsoln))
  2488.     (t lsoln))
  2489. )
  2490.  
  2491. ;;; This deals with negative, fractional, and symbolic powers of "x".
  2492. (defun math-solve-poly-funny-powers (sub-rhs)    ; uses "t1", "t2"
  2493.   (setq t1 lhs)
  2494.   (let ((pp math-poly-neg-powers)
  2495.     fac)
  2496.     (while pp
  2497.       (setq fac (math-pow (car pp) (or math-poly-mult-powers 1))
  2498.         t1 (math-mul t1 fac)
  2499.         rhs (math-mul rhs fac)
  2500.         pp (cdr pp))))
  2501.   (if sub-rhs (setq t1 (math-sub t1 rhs)))
  2502.   (let ((math-poly-neg-powers nil))
  2503.     (setq t2 (math-mul (or math-poly-mult-powers 1)
  2504.                (let ((calc-prefer-frac t))
  2505.              (math-div 1 math-poly-frac-powers)))
  2506.       t1 (math-is-polynomial (math-simplify (calcFunc-expand t1)) b 50)))
  2507. )
  2508.  
  2509. ;;; This converts "a x^8 + b x^5 + c x^2" to "(a (x^3)^2 + b (x^3) + c) * x^2".
  2510. (defun math-solve-crunch-poly (max-degree)   ; uses "t1", "t3"
  2511.   (let ((count 0))
  2512.     (while (and t1 (Math-zerop (car t1)))
  2513.       (setq t1 (cdr t1)
  2514.         count (1+ count)))
  2515.     (and t1
  2516.      (let* ((degree (1- (length t1)))
  2517.         (scale degree))
  2518.        (while (and (> scale 1) (= (car t3) 1))
  2519.          (and (= (% degree scale) 0)
  2520.           (let ((p t1)
  2521.             (n 0)
  2522.             (new-t1 nil)
  2523.             (okay t))
  2524.             (while (and p okay)
  2525.               (if (= (% n scale) 0)
  2526.               (setq new-t1 (nconc new-t1 (list (car p))))
  2527.             (or (Math-zerop (car p))
  2528.                 (setq okay nil)))
  2529.               (setq p (cdr p)
  2530.                 n (1+ n)))
  2531.             (if okay
  2532.             (setq t3 (cons scale (cdr t3))
  2533.                   t1 new-t1))))
  2534.          (setq scale (1- scale)))
  2535.        (setq t3 (list (math-mul (car t3) t2) (math-mul count t2)))
  2536.        (<= (1- (length t1)) max-degree))))
  2537. )
  2538.  
  2539. (defun calcFunc-poly (expr var &optional degree)
  2540.   (if degree
  2541.       (or (natnump degree) (math-reject-arg degree 'fixnatnump))
  2542.     (setq degree 50))
  2543.   (let ((p (math-is-polynomial expr var degree 'gen)))
  2544.     (if p
  2545.     (if (equal p '(0))
  2546.         (list 'vec)
  2547.       (cons 'vec p))
  2548.       (math-reject-arg expr "Expected a polynomial")))
  2549. )
  2550.  
  2551. (defun calcFunc-gpoly (expr var &optional degree)
  2552.   (if degree
  2553.       (or (natnump degree) (math-reject-arg degree 'fixnatnump))
  2554.     (setq degree 50))
  2555.   (let* ((math-poly-base-variable var)
  2556.      (d (math-decompose-poly expr var degree nil)))
  2557.     (if d
  2558.     (cons 'vec d)
  2559.       (math-reject-arg expr "Expected a polynomial")))
  2560. )
  2561.  
  2562. (defun math-decompose-poly (lhs solve-var degree sub-rhs)
  2563.   (let ((rhs (or sub-rhs 1))
  2564.     t1 t2 t3)
  2565.     (setq t2 (math-polynomial-base
  2566.           lhs
  2567.           (function
  2568.            (lambda (b)
  2569.          (let ((math-poly-neg-powers '(1))
  2570.                (math-poly-mult-powers nil)
  2571.                (math-poly-frac-powers 1)
  2572.                (math-poly-exp-base t))
  2573.            (and (not (equal b lhs))
  2574.             (or (not (memq (car-safe b) '(+ -))) sub-rhs)
  2575.             (setq t3 '(1 0) t2 1
  2576.                   t1 (math-is-polynomial lhs b 50))
  2577.             (if (and (equal math-poly-neg-powers '(1))
  2578.                  (memq math-poly-mult-powers '(nil 1))
  2579.                  (eq math-poly-frac-powers 1)
  2580.                  sub-rhs)
  2581.                 (setq t1 (cons (math-sub (car t1) rhs)
  2582.                        (cdr t1)))
  2583.               (math-solve-poly-funny-powers sub-rhs))
  2584.             (math-solve-crunch-poly degree)
  2585.             (or (math-expr-contains b solve-var)
  2586.                 (math-expr-contains (car t3) solve-var))))))))
  2587.     (if t2
  2588.     (list (math-pow t2 (car t3))
  2589.           (cons 'vec t1)
  2590.           (if sub-rhs
  2591.           (math-pow t2 (nth 1 t3))
  2592.         (math-div (math-pow t2 (nth 1 t3)) rhs)))))
  2593. )
  2594.  
  2595. (defun math-solve-linear (var sign b a)
  2596.   (math-try-solve-for var
  2597.               (math-div (math-neg b) a)
  2598.               (math-solve-sign sign a)
  2599.               t)
  2600. )
  2601.  
  2602. (defun math-solve-quadratic (var c b a)
  2603.   (math-try-solve-for
  2604.    var
  2605.    (if (math-looks-evenp b)
  2606.        (let ((halfb (math-div b 2)))
  2607.      (math-div
  2608.       (math-add
  2609.        (math-neg halfb)
  2610.        (math-solve-get-sign
  2611.         (math-normalize
  2612.          (list 'calcFunc-sqrt
  2613.            (math-add (math-sqr halfb)
  2614.                  (math-mul (math-neg c) a))))))
  2615.       a))
  2616.      (math-div
  2617.       (math-add
  2618.        (math-neg b)
  2619.        (math-solve-get-sign
  2620.     (math-normalize
  2621.      (list 'calcFunc-sqrt
  2622.            (math-add (math-sqr b)
  2623.              (math-mul 4 (math-mul (math-neg c) a)))))))
  2624.       (math-mul 2 a)))
  2625.    nil t)
  2626. )
  2627.  
  2628. (defun math-solve-cubic (var d c b a)
  2629.   (let* ((p (math-div b a))
  2630.      (q (math-div c a))
  2631.      (r (math-div d a))
  2632.      (psqr (math-sqr p))
  2633.      (aa (math-sub q (math-div psqr 3)))
  2634.      (bb (math-add r
  2635.                (math-div (math-sub (math-mul 2 (math-mul psqr p))
  2636.                        (math-mul 9 (math-mul p q)))
  2637.                  27)))
  2638.      m)
  2639.     (if (Math-zerop aa)
  2640.     (math-try-solve-for (math-pow (math-add var (math-div p 3)) 3)
  2641.                 (math-neg bb) nil t)
  2642.       (if (Math-zerop bb)
  2643.       (math-try-solve-for
  2644.        (math-mul (math-add var (math-div p 3))
  2645.              (math-add (math-sqr (math-add var (math-div p 3)))
  2646.                    aa))
  2647.        0 nil t)
  2648.     (setq m (math-mul 2 (list 'calcFunc-sqrt (math-div aa -3))))
  2649.     (math-try-solve-for
  2650.      var
  2651.      (math-sub
  2652.       (math-normalize
  2653.        (math-mul
  2654.         m
  2655.         (list 'calcFunc-cos
  2656.           (math-div
  2657.            (math-sub (list 'calcFunc-arccos
  2658.                    (math-div (math-mul 3 bb)
  2659.                          (math-mul aa m)))
  2660.                  (math-mul 2
  2661.                        (math-mul
  2662.                     (math-add 1 (math-solve-get-int
  2663.                              1 3))
  2664.                     (math-half-circle
  2665.                      calc-symbolic-mode))))
  2666.            3))))
  2667.       (math-div p 3))
  2668.      nil t))))
  2669. )
  2670.  
  2671. (defun math-solve-quartic (var d c b a aa)
  2672.   (setq a (math-div a aa))
  2673.   (setq b (math-div b aa))
  2674.   (setq c (math-div c aa))
  2675.   (setq d (math-div d aa))
  2676.   (math-try-solve-for
  2677.    var
  2678.    (let* ((asqr (math-sqr a))
  2679.       (asqr4 (math-div asqr 4))
  2680.       (y (let ((solve-full nil)
  2681.            calc-next-why)
  2682.            (math-solve-cubic solve-var
  2683.                  (math-sub (math-sub
  2684.                         (math-mul 4 (math-mul b d))
  2685.                         (math-mul asqr d))
  2686.                        (math-sqr c))
  2687.                  (math-sub (math-mul a c)
  2688.                        (math-mul 4 d))
  2689.                  (math-neg b)
  2690.                  1)))
  2691.       (rsqr (math-add (math-sub asqr4 b) y))
  2692.       (r (list 'calcFunc-sqrt rsqr))
  2693.       (sign1 (math-solve-get-sign 1))
  2694.       (de (list 'calcFunc-sqrt
  2695.             (math-add
  2696.              (math-sub (math-mul 3 asqr4)
  2697.                    (math-mul 2 b))
  2698.              (if (Math-zerop rsqr)
  2699.              (math-mul
  2700.               2
  2701.               (math-mul sign1
  2702.                     (list 'calcFunc-sqrt
  2703.                       (math-sub (math-sqr y)
  2704.                             (math-mul 4 d)))))
  2705.                (math-sub
  2706.             (math-mul sign1
  2707.                   (math-div
  2708.                    (math-sub (math-sub
  2709.                           (math-mul 4 (math-mul a b))
  2710.                           (math-mul 8 c))
  2711.                          (math-mul asqr a))
  2712.                    (math-mul 4 r)))
  2713.             rsqr))))))
  2714.      (math-normalize
  2715.       (math-sub (math-add (math-mul sign1 (math-div r 2))
  2716.               (math-solve-get-sign (math-div de 2)))
  2717.         (math-div a 4))))
  2718.    nil t)
  2719. )
  2720.  
  2721. (defun math-poly-all-roots (var p &optional math-factoring)
  2722.   (catch 'ouch
  2723.     (let* ((math-symbolic-solve calc-symbolic-mode)
  2724.        (roots nil)
  2725.        (deg (1- (length p)))
  2726.        (orig-p (reverse p))
  2727.        (math-int-coefs nil)
  2728.        (math-int-scale nil)
  2729.        (math-double-roots nil)
  2730.        (math-int-factors nil)
  2731.        (math-int-threshold nil)
  2732.        (pp p))
  2733.       ;; If rational coefficients, look for exact rational factors.
  2734.       (while (and pp (Math-ratp (car pp)))
  2735.     (setq pp (cdr pp)))
  2736.       (if pp
  2737.       (if (or math-factoring math-symbolic-solve)
  2738.           (throw 'ouch nil))
  2739.     (let ((lead (car orig-p))
  2740.           (calc-prefer-frac t)
  2741.           (scale (apply 'math-lcm-denoms p)))
  2742.       (setq math-int-scale (math-abs (math-mul scale lead))
  2743.         math-int-threshold (math-div '(float 5 -2) math-int-scale)
  2744.         math-int-coefs (cdr (math-div (cons 'vec orig-p) lead)))))
  2745.       (if (> deg 4)
  2746.       (let ((calc-prefer-frac nil)
  2747.         (calc-symbolic-mode nil)
  2748.         (pp p)
  2749.         (def-p (copy-sequence orig-p)))
  2750.         (while pp
  2751.           (if (Math-numberp (car pp))
  2752.           (setq pp (cdr pp))
  2753.         (throw 'ouch nil)))
  2754.         (while (> deg (if math-symbolic-solve 2 4))
  2755.           (let* ((x (math-poly-any-root def-p '(float 0 0) nil))
  2756.              b c pp)
  2757.         (if (and (eq (car-safe x) 'cplx)
  2758.              (math-nearly-zerop (nth 2 x) (nth 1 x)))
  2759.             (setq x (calcFunc-re x)))
  2760.         (or math-factoring
  2761.             (setq roots (cons x roots)))
  2762.         (or (math-numberp x)
  2763.             (setq x (math-evaluate-expr x)))
  2764.         (setq pp def-p
  2765.               b (car def-p))
  2766.         (while (setq pp (cdr pp))
  2767.           (setq c (car pp))
  2768.           (setcar pp b)
  2769.           (setq b (math-add (math-mul x b) c)))
  2770.         (setq def-p (cdr def-p)
  2771.               deg (1- deg))))
  2772.         (setq p (reverse def-p))))
  2773.       (if (> deg 1)
  2774.       (let ((solve-var '(var DUMMY var-DUMMY))
  2775.         (math-solve-sign nil)
  2776.         (math-solve-ranges nil)
  2777.         (solve-full 'all))
  2778.         (if (= (length p) (length math-int-coefs))
  2779.         (setq p (reverse math-int-coefs)))
  2780.         (setq roots (append (cdr (apply (cond ((= deg 2)
  2781.                            'math-solve-quadratic)
  2782.                           ((= deg 3)
  2783.                            'math-solve-cubic)
  2784.                           (t
  2785.                            'math-solve-quartic))
  2786.                         solve-var p))
  2787.                 roots)))
  2788.     (if (> deg 0)
  2789.         (setq roots (cons (math-div (math-neg (car p)) (nth 1 p))
  2790.                   roots))))
  2791.       (if math-factoring
  2792.       (progn
  2793.         (while roots
  2794.           (math-poly-integer-root (car roots))
  2795.           (setq roots (cdr roots)))
  2796.         (list math-int-factors (nreverse math-int-coefs) math-int-scale))
  2797.     (let ((vec nil) res)
  2798.       (while roots
  2799.         (let ((root (car roots))
  2800.           (solve-full (and solve-full 'all)))
  2801.           (if (math-floatp root)
  2802.           (setq root (math-poly-any-root orig-p root t)))
  2803.           (setq vec (append vec
  2804.                 (cdr (or (math-try-solve-for var root nil t)
  2805.                      (throw 'ouch nil))))))
  2806.         (setq roots (cdr roots)))
  2807.       (setq vec (cons 'vec (nreverse vec)))
  2808.       (if math-symbolic-solve
  2809.           (setq vec (math-normalize vec)))
  2810.       (if (eq solve-full t)
  2811.           (list 'calcFunc-subscr
  2812.             vec
  2813.             (math-solve-get-int 1 (1- (length orig-p)) 1))
  2814.         vec)))))
  2815. )
  2816. (setq math-symbolic-solve nil)
  2817.  
  2818. (defun math-lcm-denoms (&rest fracs)
  2819.   (let ((den 1))
  2820.     (while fracs
  2821.       (if (eq (car-safe (car fracs)) 'frac)
  2822.       (setq den (calcFunc-lcm den (nth 2 (car fracs)))))
  2823.       (setq fracs (cdr fracs)))
  2824.     den)
  2825. )
  2826.  
  2827. (defun math-poly-any-root (p x polish)    ; p is a reverse poly coeff list
  2828.   (let* ((newt (if (math-zerop x)
  2829.            (math-poly-newton-root
  2830.             p '(cplx (float 123 -6) (float 1 -4)) 4)
  2831.          (math-poly-newton-root p x 4)))
  2832.      (res (if (math-zerop (cdr newt))
  2833.           (car newt)
  2834.         (if (and (math-lessp (cdr newt) '(float 1 -3)) (not polish))
  2835.             (setq newt (math-poly-newton-root p (car newt) 30)))
  2836.         (if (math-zerop (cdr newt))
  2837.             (car newt)
  2838.           (math-poly-laguerre-root p x polish)))))
  2839.     (and math-symbolic-solve (math-floatp res)
  2840.      (throw 'ouch nil))
  2841.     res)
  2842. )
  2843.  
  2844. (defun math-poly-newton-root (p x iters)
  2845.   (let* ((calc-prefer-frac nil)
  2846.      (calc-symbolic-mode nil)
  2847.      (try-integer math-int-coefs)
  2848.      (dx x) b d)
  2849.     (while (and (> (setq iters (1- iters)) 0)
  2850.         (let ((pp p))
  2851.           (math-working "newton" x)
  2852.           (setq b (car p)
  2853.             d 0)
  2854.           (while (setq pp (cdr pp))
  2855.             (setq d (math-add (math-mul x d) b)
  2856.               b (math-add (math-mul x b) (car pp))))
  2857.           (not (math-zerop d)))
  2858.         (progn
  2859.           (setq dx (math-div b d)
  2860.             x (math-sub x dx))
  2861.           (if try-integer
  2862.               (let ((adx (math-abs-approx dx)))
  2863.             (and (math-lessp adx math-int-threshold)
  2864.                  (let ((iroot (math-poly-integer-root x)))
  2865.                    (if iroot
  2866.                    (setq x iroot dx 0)
  2867.                  (setq try-integer nil))))))
  2868.           (or (not (or (eq dx 0)
  2869.                    (math-nearly-zerop dx (math-abs-approx x))))
  2870.               (progn (setq dx 0) nil)))))
  2871.     (cons x (if (math-zerop x)
  2872.         1 (math-div (math-abs-approx dx) (math-abs-approx x)))))
  2873. )
  2874.  
  2875. (defun math-poly-integer-root (x)
  2876.   (and (math-lessp (calcFunc-xpon (math-abs-approx x)) calc-internal-prec)
  2877.        math-int-coefs
  2878.        (let* ((calc-prefer-frac t)
  2879.           (xre (calcFunc-re x))
  2880.           (xim (calcFunc-im x))
  2881.           (xresq (math-sqr xre))
  2882.           (ximsq (math-sqr xim)))
  2883.      (if (math-lessp ximsq (calcFunc-scf xresq -1))
  2884.          ;; Look for linear factor
  2885.          (let* ((rnd (math-div (math-round (math-mul xre math-int-scale))
  2886.                    math-int-scale))
  2887.             (icp math-int-coefs)
  2888.             (rem (car icp))
  2889.             (newcoef nil))
  2890.            (while (setq icp (cdr icp))
  2891.          (setq newcoef (cons rem newcoef)
  2892.                rem (math-add (car icp)
  2893.                      (math-mul rem rnd))))
  2894.            (and (math-zerop rem)
  2895.             (progn
  2896.               (setq math-int-coefs (nreverse newcoef)
  2897.                 math-int-factors (cons (list (math-neg rnd))
  2898.                            math-int-factors))
  2899.               rnd)))
  2900.        ;; Look for irreducible quadratic factor
  2901.        (let* ((rnd1 (math-div (math-round
  2902.                    (math-mul xre (math-mul -2 math-int-scale)))
  2903.                   math-int-scale))
  2904.           (sqscale (math-sqr math-int-scale))
  2905.           (rnd0 (math-div (math-round (math-mul (math-add xresq ximsq)
  2906.                             sqscale))
  2907.                   sqscale))
  2908.           (rem1 (car math-int-coefs))
  2909.           (icp (cdr math-int-coefs))
  2910.           (rem0 (car icp))
  2911.           (newcoef nil)
  2912.           (found (assoc (list rnd0 rnd1 (math-posp xim))
  2913.                 math-double-roots))
  2914.           this)
  2915.          (if found
  2916.          (setq math-double-roots (delq found math-double-roots)
  2917.                rem0 0 rem1 0)
  2918.            (while (setq icp (cdr icp))
  2919.          (setq this rem1
  2920.                newcoef (cons rem1 newcoef)
  2921.                rem1 (math-sub rem0 (math-mul this rnd1))
  2922.                rem0 (math-sub (car icp) (math-mul this rnd0)))))
  2923.          (and (math-zerop rem0)
  2924.           (math-zerop rem1)
  2925.           (let ((aa (math-div rnd1 -2)))
  2926.             (or found (setq math-int-coefs (reverse newcoef)
  2927.                     math-double-roots (cons (list
  2928.                                  (list
  2929.                                   rnd0 rnd1
  2930.                                   (math-negp xim)))
  2931.                                 math-double-roots)
  2932.                     math-int-factors (cons (cons rnd0 rnd1)
  2933.                                math-int-factors)))
  2934.             (math-add aa
  2935.                   (let ((calc-symbolic-mode math-symbolic-solve))
  2936.                 (math-mul (math-sqrt (math-sub (math-sqr aa)
  2937.                                    rnd0))
  2938.                       (if (math-negp xim) -1 1))))))))))
  2939. )
  2940. (setq math-int-coefs nil)
  2941.  
  2942. ;;; The following routine is from Numerical Recipes, section 9.5.
  2943. (defun math-poly-laguerre-root (p x polish)
  2944.   (let* ((calc-prefer-frac nil)
  2945.      (calc-symbolic-mode nil)
  2946.      (iters 0)
  2947.      (m (1- (length p)))
  2948.      (try-newt (not polish))
  2949.      (tried-newt nil)
  2950.      b d f x1 dx dxold)
  2951.     (while
  2952.     (and (or (< (setq iters (1+ iters)) 50)
  2953.          (math-reject-arg x "*Laguerre's method failed to converge"))
  2954.          (let ((err (math-abs-approx (car p)))
  2955.            (abx (math-abs-approx x))
  2956.            (pp p))
  2957.            (setq b (car p)
  2958.              d 0 f 0)
  2959.            (while (setq pp (cdr pp))
  2960.          (setq f (math-add (math-mul x f) d)
  2961.                d (math-add (math-mul x d) b)
  2962.                b (math-add (math-mul x b) (car pp))
  2963.                err (math-add (math-abs-approx b) (math-mul abx err))))
  2964.            (math-lessp (calcFunc-scf err (- -2 calc-internal-prec))
  2965.                (math-abs-approx b)))
  2966.          (or (not (math-zerop d))
  2967.          (not (math-zerop f))
  2968.          (progn
  2969.            (setq x (math-pow (math-neg b) (list 'frac 1 m)))
  2970.            nil))
  2971.          (let* ((g (math-div d b))
  2972.             (g2 (math-sqr g))
  2973.             (h (math-sub g2 (math-mul 2 (math-div f b))))
  2974.             (sq (math-sqrt
  2975.              (math-mul (1- m) (math-sub (math-mul m h) g2))))
  2976.             (gp (math-add g sq))
  2977.             (gm (math-sub g sq)))
  2978.            (if (math-lessp (calcFunc-abssqr gp) (calcFunc-abssqr gm))
  2979.            (setq gp gm))
  2980.            (setq dx (math-div m gp)
  2981.              x1 (math-sub x dx))
  2982.            (if (and try-newt
  2983.             (math-lessp (math-abs-approx dx)
  2984.                     (calcFunc-scf (math-abs-approx x) -3)))
  2985.            (let ((newt (math-poly-newton-root p x1 7)))
  2986.              (setq tried-newt t
  2987.                try-newt nil)
  2988.              (if (math-zerop (cdr newt))
  2989.              (setq x (car newt) x1 x)
  2990.                (if (math-lessp (cdr newt) '(float 1 -6))
  2991.                (let ((newt2 (math-poly-newton-root
  2992.                      p (car newt) 20)))
  2993.                  (if (math-zerop (cdr newt2))
  2994.                  (setq x (car newt2) x1 x)
  2995.                    (setq x (car newt))))))))
  2996.            (not (or (eq x x1)
  2997.             (math-nearly-equal x x1))))
  2998.          (let ((cdx (math-abs-approx dx)))
  2999.            (setq x x1
  3000.              tried-newt nil)
  3001.            (prog1
  3002.            (or (<= iters 6)
  3003.                (math-lessp cdx dxold)
  3004.                (progn
  3005.              (if polish
  3006.                  (let ((digs (calcFunc-xpon
  3007.                       (math-div (math-abs-approx x) cdx))))
  3008.                    (calc-record-why
  3009.                 "*Could not attain full precision")
  3010.                    (if (natnump digs)
  3011.                    (let ((calc-internal-prec (max 3 digs)))
  3012.                      (setq x (math-normalize x))))))
  3013.              nil))
  3014.          (setq dxold cdx)))
  3015.          (or polish
  3016.          (math-lessp (calcFunc-scf (math-abs-approx x)
  3017.                        (- calc-internal-prec))
  3018.                  dxold))))
  3019.     (or (and (math-floatp x)
  3020.          (math-poly-integer-root x))
  3021.     x))
  3022. )
  3023.  
  3024. (defun math-solve-above-dummy (x)
  3025.   (and (not (Math-primp x))
  3026.        (if (and (equal (nth 1 x) '(var SOLVEDUM SOLVEDUM))
  3027.         (= (length x) 2))
  3028.        x
  3029.      (let ((res nil))
  3030.        (while (and (setq x (cdr x))
  3031.                (not (setq res (math-solve-above-dummy (car x))))))
  3032.        res)))
  3033. )
  3034.  
  3035. (defun math-solve-find-root-term (x neg)    ; sets "t2", "t3"
  3036.   (if (math-solve-find-root-in-prod x)
  3037.       (setq t3 neg
  3038.         t1 x)
  3039.     (and (memq (car-safe x) '(+ -))
  3040.      (or (math-solve-find-root-term (nth 1 x) neg)
  3041.          (math-solve-find-root-term (nth 2 x)
  3042.                     (if (eq (car x) '-) (not neg) neg)))))
  3043. )
  3044.  
  3045. (defun math-solve-find-root-in-prod (x)
  3046.   (and (consp x)
  3047.        (math-expr-contains x solve-var)
  3048.        (or (and (eq (car x) 'calcFunc-sqrt)
  3049.         (setq t2 2))
  3050.        (and (eq (car x) '^)
  3051.         (or (and (memq (math-quarter-integer (nth 2 x)) '(1 2 3))
  3052.              (setq t2 2))
  3053.             (and (eq (car-safe (nth 2 x)) 'frac)
  3054.              (eq (nth 2 (nth 2 x)) 3)
  3055.              (setq t2 3))))
  3056.        (and (memq (car x) '(* /))
  3057.         (or (and (not (math-expr-contains (nth 1 x) solve-var))
  3058.              (math-solve-find-root-in-prod (nth 2 x)))
  3059.             (and (not (math-expr-contains (nth 2 x) solve-var))
  3060.              (math-solve-find-root-in-prod (nth 1 x)))))))
  3061. )
  3062.  
  3063.  
  3064. (defun math-solve-system (exprs solve-vars solve-full)
  3065.   (setq exprs (mapcar 'list (if (Math-vectorp exprs)
  3066.                 (cdr exprs)
  3067.                   (list exprs)))
  3068.     solve-vars (if (Math-vectorp solve-vars)
  3069.                (cdr solve-vars)
  3070.              (list solve-vars)))
  3071.   (or (let ((math-solve-simplifying nil))
  3072.     (math-solve-system-rec exprs solve-vars nil))
  3073.       (let ((math-solve-simplifying t))
  3074.     (math-solve-system-rec exprs solve-vars nil)))
  3075. )
  3076.  
  3077. ;;; The following backtracking solver works by choosing a variable
  3078. ;;; and equation, and trying to solve the equation for the variable.
  3079. ;;; If it succeeds it calls itself recursively with that variable and
  3080. ;;; equation removed from their respective lists, and with the solution
  3081. ;;; added to solns as well as being substituted into all existing
  3082. ;;; equations.  The algorithm terminates when any solution path
  3083. ;;; manages to remove all the variables from var-list.
  3084.  
  3085. ;;; To support calcFunc-roots, entries in eqn-list and solns are
  3086. ;;; actually lists of equations.
  3087.  
  3088. (defun math-solve-system-rec (eqn-list var-list solns)
  3089.   (if var-list
  3090.       (let ((v var-list)
  3091.         (res nil))
  3092.  
  3093.     ;; Try each variable in turn.
  3094.     (while
  3095.         (and
  3096.          v
  3097.          (let* ((vv (car v))
  3098.             (e eqn-list)
  3099.             (elim (eq (car-safe vv) 'calcFunc-elim)))
  3100.            (if elim
  3101.            (setq vv (nth 1 vv)))
  3102.  
  3103.            ;; Try each equation in turn.
  3104.            (while
  3105.            (and
  3106.             e
  3107.             (let ((e2 (car e))
  3108.               (eprev nil)
  3109.               res2)
  3110.               (setq res nil)
  3111.  
  3112.               ;; Try to solve for vv the list of equations e2.
  3113.               (while (and e2
  3114.                   (setq res2 (or (and (eq (car e2) eprev)
  3115.                               res2)
  3116.                          (math-solve-for (car e2) 0 vv
  3117.                                  solve-full))))
  3118.             (setq eprev (car e2)
  3119.                   res (cons (if (eq solve-full 'all)
  3120.                         (cdr res2)
  3121.                       (list res2))
  3122.                     res)
  3123.                   e2 (cdr e2)))
  3124.               (if e2
  3125.               (setq res nil)
  3126.  
  3127.             ;; Found a solution.  Now try other variables.
  3128.             (setq res (nreverse res)
  3129.                   res (math-solve-system-rec
  3130.                    (mapcar
  3131.                     'math-solve-system-subst
  3132.                     (delq (car e)
  3133.                       (copy-sequence eqn-list)))
  3134.                    (delq (car v) (copy-sequence var-list))
  3135.                    (let ((math-solve-simplifying nil)
  3136.                      (s (mapcar
  3137.                          (function
  3138.                           (lambda (x)
  3139.                         (cons
  3140.                          (car x)
  3141.                          (math-solve-system-subst
  3142.                           (cdr x)))))
  3143.                          solns)))
  3144.                      (if elim
  3145.                      s
  3146.                        (cons (cons vv (apply 'append res))
  3147.                          s)))))
  3148.             (not res))))
  3149.          (setq e (cdr e)))
  3150.            (not res)))
  3151.       (setq v (cdr v)))
  3152.     res)
  3153.  
  3154.     ;; Eliminated all variables, so now put solution into the proper format.
  3155.     (setq solns (sort solns
  3156.               (function
  3157.                (lambda (x y)
  3158.              (not (memq (car x) (memq (car y) solve-vars)))))))
  3159.     (if (eq solve-full 'all)
  3160.     (math-transpose
  3161.      (math-normalize
  3162.       (cons 'vec
  3163.         (if solns
  3164.             (mapcar (function (lambda (x) (cons 'vec (cdr x)))) solns)
  3165.           (mapcar (function (lambda (x) (cons 'vec x))) eqn-list)))))
  3166.       (math-normalize
  3167.        (cons 'vec 
  3168.          (if solns
  3169.          (mapcar (function (lambda (x) (cons 'calcFunc-eq x))) solns)
  3170.            (mapcar 'car eqn-list))))))
  3171. )
  3172.  
  3173. (defun math-solve-system-subst (x)    ; uses "res" and "v"
  3174.   (let ((accum nil)
  3175.     (res2 res))
  3176.     (while x
  3177.       (setq accum (nconc accum
  3178.              (mapcar (function
  3179.                   (lambda (r)
  3180.                     (if math-solve-simplifying
  3181.                     (math-simplify
  3182.                      (math-expr-subst (car x) vv r))
  3183.                       (math-expr-subst (car x) vv r))))
  3184.                  (car res2)))
  3185.         x (cdr x)
  3186.         res2 (cdr res2)))
  3187.     accum)
  3188. )
  3189.  
  3190.  
  3191. (defun math-get-from-counter (name)
  3192.   (let ((ctr (assq name calc-command-flags)))
  3193.     (if ctr
  3194.     (setcdr ctr (1+ (cdr ctr)))
  3195.       (setq ctr (cons name 1)
  3196.         calc-command-flags (cons ctr calc-command-flags)))
  3197.     (cdr ctr))
  3198. )
  3199.  
  3200. (defun math-solve-get-sign (val)
  3201.   (setq val (math-simplify val))
  3202.   (if (and (eq (car-safe val) '*)
  3203.        (Math-numberp (nth 1 val)))
  3204.       (list '* (nth 1 val) (math-solve-get-sign (nth 2 val)))
  3205.     (and (eq (car-safe val) 'calcFunc-sqrt)
  3206.      (eq (car-safe (nth 1 val)) '^)
  3207.      (setq val (math-normalize (list '^
  3208.                      (nth 1 (nth 1 val))
  3209.                      (math-div (nth 2 (nth 1 val)) 2)))))
  3210.     (if solve-full
  3211.     (if (and (calc-var-value 'var-GenCount)
  3212.          (Math-natnump var-GenCount)
  3213.          (not (eq solve-full 'all)))
  3214.         (prog1
  3215.         (math-mul (list 'calcFunc-as var-GenCount) val)
  3216.           (setq var-GenCount (math-add var-GenCount 1))
  3217.           (calc-refresh-evaltos 'var-GenCount))
  3218.       (let* ((var (concat "s" (math-get-from-counter 'solve-sign)))
  3219.          (var2 (list 'var (intern var) (intern (concat "var-" var)))))
  3220.         (if (eq solve-full 'all)
  3221.         (setq math-solve-ranges (cons (list var2 1 -1)
  3222.                           math-solve-ranges)))
  3223.         (math-mul var2 val)))
  3224.       (calc-record-why "*Choosing positive solution")
  3225.       val))
  3226. )
  3227.  
  3228. (defun math-solve-get-int (val &optional range first)
  3229.   (if solve-full
  3230.       (if (and (calc-var-value 'var-GenCount)
  3231.            (Math-natnump var-GenCount)
  3232.            (not (eq solve-full 'all)))
  3233.       (prog1
  3234.           (math-mul val (list 'calcFunc-an var-GenCount))
  3235.         (setq var-GenCount (math-add var-GenCount 1))
  3236.         (calc-refresh-evaltos 'var-GenCount))
  3237.     (let* ((var (concat "n" (math-get-from-counter 'solve-int)))
  3238.            (var2 (list 'var (intern var) (intern (concat "var-" var)))))
  3239.       (if (and range (eq solve-full 'all))
  3240.           (setq math-solve-ranges (cons (cons var2
  3241.                           (cdr (calcFunc-index
  3242.                             range (or first 0))))
  3243.                         math-solve-ranges)))
  3244.       (math-mul val var2)))
  3245.     (calc-record-why "*Choosing 0 for arbitrary integer in solution")
  3246.     0)
  3247. )
  3248.  
  3249. (defun math-solve-sign (sign expr)
  3250.   (and sign
  3251.        (let ((s1 (math-possible-signs expr)))
  3252.      (cond ((memq s1 '(4 6))
  3253.         sign)
  3254.            ((memq s1 '(1 3))
  3255.         (- sign)))))
  3256. )
  3257.  
  3258. (defun math-looks-evenp (expr)
  3259.   (if (Math-integerp expr)
  3260.       (math-evenp expr)
  3261.     (if (memq (car expr) '(* /))
  3262.     (math-looks-evenp (nth 1 expr))))
  3263. )
  3264.  
  3265. (defun math-solve-for (lhs rhs solve-var solve-full &optional sign)
  3266.   (if (math-expr-contains rhs solve-var)
  3267.       (math-solve-for (math-sub lhs rhs) 0 solve-var solve-full)
  3268.     (and (math-expr-contains lhs solve-var)
  3269.      (math-with-extra-prec 1
  3270.        (let* ((math-poly-base-variable solve-var)
  3271.           (res (math-try-solve-for lhs rhs sign)))
  3272.          (if (and (eq solve-full 'all)
  3273.               (math-known-realp solve-var))
  3274.          (let ((old-len (length res))
  3275.                new-len)
  3276.            (setq res (delq nil
  3277.                    (mapcar (function
  3278.                         (lambda (x)
  3279.                           (and (not (memq (car-safe x)
  3280.                                   '(cplx polar)))
  3281.                            x)))
  3282.                        res))
  3283.              new-len (length res))
  3284.            (if (< new-len old-len)
  3285.                (calc-record-why (if (= new-len 1)
  3286.                         "*All solutions were complex"
  3287.                       (format
  3288.                        "*Omitted %d complex solutions"
  3289.                        (- old-len new-len)))))))
  3290.          res))))
  3291. )
  3292.  
  3293. (defun math-solve-eqn (expr var full)
  3294.   (if (memq (car-safe expr) '(calcFunc-neq calcFunc-lt calcFunc-gt
  3295.                        calcFunc-leq calcFunc-geq))
  3296.       (let ((res (math-solve-for (cons '- (cdr expr))
  3297.                  0 var full
  3298.                  (if (eq (car expr) 'calcFunc-neq) nil 1))))
  3299.     (and res
  3300.          (if (eq math-solve-sign 1)
  3301.          (list (car expr) var res)
  3302.            (if (eq math-solve-sign -1)
  3303.            (list (car expr) res var)
  3304.          (or (eq (car expr) 'calcFunc-neq)
  3305.              (calc-record-why
  3306.               "*Can't determine direction of inequality"))
  3307.          (and (memq (car expr) '(calcFunc-neq calcFunc-lt calcFunc-gt))
  3308.               (list 'calcFunc-neq var res))))))
  3309.     (let ((res (math-solve-for expr 0 var full)))
  3310.       (and res
  3311.        (list 'calcFunc-eq var res))))
  3312. )
  3313.  
  3314. (defun math-reject-solution (expr var func)
  3315.   (if (math-expr-contains expr var)
  3316.       (or (equal (car calc-next-why) '(* "Unable to find a symbolic solution"))
  3317.       (calc-record-why "*Unable to find a solution")))
  3318.   (list func expr var)
  3319. )
  3320.  
  3321. (defun calcFunc-solve (expr var)
  3322.   (or (if (or (Math-vectorp expr) (Math-vectorp var))
  3323.       (math-solve-system expr var nil)
  3324.     (math-solve-eqn expr var nil))
  3325.       (math-reject-solution expr var 'calcFunc-solve))
  3326. )
  3327.  
  3328. (defun calcFunc-fsolve (expr var)
  3329.   (or (if (or (Math-vectorp expr) (Math-vectorp var))
  3330.       (math-solve-system expr var t)
  3331.     (math-solve-eqn expr var t))
  3332.       (math-reject-solution expr var 'calcFunc-fsolve))
  3333. )
  3334.  
  3335. (defun calcFunc-roots (expr var)
  3336.   (let ((math-solve-ranges nil))
  3337.     (or (if (or (Math-vectorp expr) (Math-vectorp var))
  3338.         (math-solve-system expr var 'all)
  3339.       (math-solve-for expr 0 var 'all))
  3340.       (math-reject-solution expr var 'calcFunc-roots)))
  3341. )
  3342.  
  3343. (defun calcFunc-finv (expr var)
  3344.   (let ((res (math-solve-for expr math-integ-var var nil)))
  3345.     (if res
  3346.     (math-normalize (math-expr-subst res math-integ-var var))
  3347.       (math-reject-solution expr var 'calcFunc-finv)))
  3348. )
  3349.  
  3350. (defun calcFunc-ffinv (expr var)
  3351.   (let ((res (math-solve-for expr math-integ-var var t)))
  3352.     (if res
  3353.     (math-normalize (math-expr-subst res math-integ-var var))
  3354.       (math-reject-solution expr var 'calcFunc-finv)))
  3355. )
  3356.  
  3357.  
  3358. (put 'calcFunc-inv 'math-inverse
  3359.      (function (lambda (x) (math-div 1 x))))
  3360. (put 'calcFunc-inv 'math-inverse-sign -1)
  3361.  
  3362. (put 'calcFunc-sqrt 'math-inverse
  3363.      (function (lambda (x) (math-sqr x))))
  3364.  
  3365. (put 'calcFunc-conj 'math-inverse
  3366.      (function (lambda (x) (list 'calcFunc-conj x))))
  3367.  
  3368. (put 'calcFunc-abs 'math-inverse
  3369.      (function (lambda (x) (math-solve-get-sign x))))
  3370.  
  3371. (put 'calcFunc-deg 'math-inverse
  3372.      (function (lambda (x) (list 'calcFunc-rad x))))
  3373. (put 'calcFunc-deg 'math-inverse-sign 1)
  3374.  
  3375. (put 'calcFunc-rad 'math-inverse
  3376.      (function (lambda (x) (list 'calcFunc-deg x))))
  3377. (put 'calcFunc-rad 'math-inverse-sign 1)
  3378.  
  3379. (put 'calcFunc-ln 'math-inverse
  3380.      (function (lambda (x) (list 'calcFunc-exp x))))
  3381. (put 'calcFunc-ln 'math-inverse-sign 1)
  3382.  
  3383. (put 'calcFunc-log10 'math-inverse
  3384.      (function (lambda (x) (list 'calcFunc-exp10 x))))
  3385. (put 'calcFunc-log10 'math-inverse-sign 1)
  3386.  
  3387. (put 'calcFunc-lnp1 'math-inverse
  3388.      (function (lambda (x) (list 'calcFunc-expm1 x))))
  3389. (put 'calcFunc-lnp1 'math-inverse-sign 1)
  3390.  
  3391. (put 'calcFunc-exp 'math-inverse
  3392.      (function (lambda (x) (math-add (math-normalize (list 'calcFunc-ln x))
  3393.                      (math-mul 2
  3394.                            (math-mul '(var pi var-pi)
  3395.                              (math-solve-get-int
  3396.                               '(var i var-i))))))))
  3397. (put 'calcFunc-exp 'math-inverse-sign 1)
  3398.  
  3399. (put 'calcFunc-expm1 'math-inverse
  3400.      (function (lambda (x) (math-add (math-normalize (list 'calcFunc-lnp1 x))
  3401.                      (math-mul 2
  3402.                            (math-mul '(var pi var-pi)
  3403.                              (math-solve-get-int
  3404.                               '(var i var-i))))))))
  3405. (put 'calcFunc-expm1 'math-inverse-sign 1)
  3406.  
  3407. (put 'calcFunc-sin 'math-inverse
  3408.      (function (lambda (x) (let ((n (math-solve-get-int 1)))
  3409.                  (math-add (math-mul (math-normalize
  3410.                           (list 'calcFunc-arcsin x))
  3411.                          (math-pow -1 n))
  3412.                        (math-mul (math-half-circle t)
  3413.                          n))))))
  3414.  
  3415. (put 'calcFunc-cos 'math-inverse
  3416.      (function (lambda (x) (math-add (math-solve-get-sign
  3417.                       (math-normalize
  3418.                        (list 'calcFunc-arccos x)))
  3419.                      (math-solve-get-int
  3420.                       (math-full-circle t))))))
  3421.  
  3422. (put 'calcFunc-tan 'math-inverse
  3423.      (function (lambda (x) (math-add (math-normalize (list 'calcFunc-arctan x))
  3424.                      (math-solve-get-int
  3425.                       (math-half-circle t))))))
  3426.  
  3427. (put 'calcFunc-arcsin 'math-inverse
  3428.      (function (lambda (x) (math-normalize (list 'calcFunc-sin x)))))
  3429.  
  3430. (put 'calcFunc-arccos 'math-inverse
  3431.      (function (lambda (x) (math-normalize (list 'calcFunc-cos x)))))
  3432.  
  3433. (put 'calcFunc-arctan 'math-inverse
  3434.      (function (lambda (x) (math-normalize (list 'calcFunc-tan x)))))
  3435.  
  3436. (put 'calcFunc-sinh 'math-inverse
  3437.      (function (lambda (x) (let ((n (math-solve-get-int 1)))
  3438.                  (math-add (math-mul (math-normalize
  3439.                           (list 'calcFunc-arcsinh x))
  3440.                          (math-pow -1 n))
  3441.                        (math-mul (math-half-circle t)
  3442.                          (math-mul
  3443.                           '(var i var-i)
  3444.                           n)))))))
  3445. (put 'calcFunc-sinh 'math-inverse-sign 1)
  3446.  
  3447. (put 'calcFunc-cosh 'math-inverse
  3448.      (function (lambda (x) (math-add (math-solve-get-sign
  3449.                       (math-normalize
  3450.                        (list 'calcFunc-arccosh x)))
  3451.                      (math-mul (math-full-circle t)
  3452.                            (math-solve-get-int
  3453.                         '(var i var-i)))))))
  3454.  
  3455. (put 'calcFunc-tanh 'math-inverse
  3456.      (function (lambda (x) (math-add (math-normalize
  3457.                       (list 'calcFunc-arctanh x))
  3458.                      (math-mul (math-half-circle t)
  3459.                            (math-solve-get-int
  3460.                         '(var i var-i)))))))
  3461. (put 'calcFunc-tanh 'math-inverse-sign 1)
  3462.  
  3463. (put 'calcFunc-arcsinh 'math-inverse
  3464.      (function (lambda (x) (math-normalize (list 'calcFunc-sinh x)))))
  3465. (put 'calcFunc-arcsinh 'math-inverse-sign 1)
  3466.  
  3467. (put 'calcFunc-arccosh 'math-inverse
  3468.      (function (lambda (x) (math-normalize (list 'calcFunc-cosh x)))))
  3469.  
  3470. (put 'calcFunc-arctanh 'math-inverse
  3471.      (function (lambda (x) (math-normalize (list 'calcFunc-tanh x)))))
  3472. (put 'calcFunc-arctanh 'math-inverse-sign 1)
  3473.  
  3474.  
  3475.  
  3476. (defun calcFunc-taylor (expr var num)
  3477.   (let ((x0 0) (v var))
  3478.     (if (memq (car-safe var) '(+ - calcFunc-eq))
  3479.     (setq x0 (if (eq (car var) '+) (math-neg (nth 2 var)) (nth 2 var))
  3480.           v (nth 1 var)))
  3481.     (or (and (eq (car-safe v) 'var)
  3482.          (math-expr-contains expr v)
  3483.          (natnump num)
  3484.          (let ((accum (math-expr-subst expr v x0))
  3485.            (var2 (if (eq (car var) 'calcFunc-eq)
  3486.                  (cons '- (cdr var))
  3487.                var))
  3488.            (n 0)
  3489.            (nfac 1)
  3490.            (fprime expr))
  3491.            (while (and (<= (setq n (1+ n)) num)
  3492.                (setq fprime (calcFunc-deriv fprime v nil t)))
  3493.          (setq fprime (math-simplify fprime)
  3494.                nfac (math-mul nfac n)
  3495.                accum (math-add accum
  3496.                        (math-div (math-mul (math-pow var2 n)
  3497.                                (math-expr-subst
  3498.                                 fprime v x0))
  3499.                          nfac))))
  3500.            (and fprime
  3501.             (math-normalize accum))))
  3502.     (list 'calcFunc-taylor expr var num)))
  3503. )
  3504.  
  3505.  
  3506.  
  3507.  
  3508.